In Vitro Evaluation of the Effects of Prebiotics GOS and FOS on Human Colonic Caco – 2 Cells

: Prebiotic oligosaccharides are widely used as human and animal feed additives for their beneficial effects on the gut microbiota. However, there are limited data to assess the direct effect of such functional foods on the transcriptome of intestinal epithelial cells. The purpose of this study is to describe the differential transcriptomes and cellular pathways of colonic cells directly exposed to galactooligosaccharides (GOS) and fructo-oligosaccharides (FOS). We have examined the differential gene expression of polarized Caco-2 cells treated with GOS or FOS and their respective mock-treated cells using mRNA sequencing (RNA-seq). A total of 89 significant differentially expressed genes were identified between GOS and mock-treated groups. For FOS treatment, a reduced number of 12 significant genes were observed to be differentially expressed relative to the control group. KEGG and Gene Ontology functional analysis revealed that genes up-regulated in the presence of GOS were involved in digestion and absorption processes, fatty acids and steroids metabolism, potential antimicrobial proteins, energy-dependent and independent transmembrane trafficking of solutes and amino acids. Using our data, we have established complementary non-prebiotic modes of action for these frequently used dietary fibers.


Introduction
Prebiotics are generally defined as substances that are selectively utilized by host micro-organisms conferring a health benefit [1]. Prebiotic oligosaccharides conform to the definition as non-digestible food ingredients that beneficially affect the host by stimulating the growth and/or activity of beneficial bacteria in the colon [2]. Such oligosaccharides are widely used as human and animal nutritional additives for their beneficial effects on the composition of the probiotic microbiota and gut health [3,4]. However, oligosaccharides that largely resist hydrolytic activities of salivary and intestinal digestive enzymes reach the colon virtually intact, thus prompting to question whether they could have a direct or 'non-prebiotic' effect on the intestinal epithelium should they reach the mucosa.
In vivo, galacto-oligosaccharides supplementation has been shown to increase iron absorption from micronutrient powder in Kenyan infants [8] or to reduce stressinduced gastrointestinal dysfunction and days of cold or flu symptoms in controlled trials of healthy university students [9]. In a healthy volunteer study, Schmidt et al. [10] described the reduction of the salivary cortisol awakening response associated with GOS administration, implying suppression of the neuroendocrine stress response. However, Liu et al. [11] suggested short-term intake of high-dose GOS and FOS prebiotics had an adverse effect on glucose metabolism despite increased Bifidobacterium in the fecal microbiota. In controlled animal studies, GOS has been associated with improved growth performance of broiler chickens [12] or during exposure to heat stress [13,14]. Dietary supplementation of GOS has also been linked with improvements in the transition to a mature intestinal microbiota in broiler chickens [12,15] and in suckling piglets [16,17]. In a study conducted with suckling male rats, in addition to changes in the intestinal microbiota, Le Dréan et al. [18] observed GOS and FOS supplementation impacted entero-endocrine cell maturation by bringing about transient increases in the density of GLP-1 cells and the production of the satiety-related peptides GLP-1 and PYY.
The direct cellular effects of prebiotics have been investigated using polarized intestinal epithelial cell models. For example, GOS exposure has been associated with improved tight junction formation and re-epithelialization post disruption, which coincided with increases in the expression of cell proliferation and cell differentiation pathways [19]. GOS was also reported to prevent deoxynivalenolinduced compromise of the integrity of Caco-2 cell monolayers and the corresponding decline in the expression of the tight junction encoding gene CLDN3 [20]. GOS is also reported to reduce the adherence of Salmonella Typhimurium to mucus-and to non-mucus-producing HT-29 cells [21].
Several studies have reported that GOS, FOS and inulin promote calcium absorption in animal and human gut studies [22][23][24]. FOS has also been reported to have non-microbiome or non-prebiotic mediated effects on immune function. For example, Fransen et al. used germ-free mice to confirm microbiota-independent changes in immune function with short-or long-chain β2→1-fructans enhancing T helper 1 cells in Peyer's patches and short-chain β2→1-fructans increasing 3 of 34 regulatory T cells and CD11b−CD103− dendritic cells in the mesenteric lymph nodes [25]. FOS-inulin pre-incubation of a chicken macrophage cell line HD11 before challenge with Salmonella Enteriditis reduced cellular uptake of the pathogen and IL-1β gene expression suggesting that inulin-enriched FOS had the ability to modulate the innate immune response [26].
Caco-2 cells are a continuous line of heterogeneous human epithelial colorectal adenocarcinoma cells that produce tight junctions, microvilli, and enzymes and transporters characteristic of enterocytes. The Caco-2 monolayer is widely used in the pharmaceutical industry as an in vitro model of the human intestinal mucosa to predict the absorption of orally administered drugs. Caco-2 cells are commonly cultured as a polarized epithelial stratum by forming confluent monolayers on insert filters that provide a physical and biochemical barrier to the passage of ions and small molecules [27]. In this study we have examined the effects of Nutrabiotic® GOS and Beneo® FOS (also referred as "GOS" and "FOS" respectively) on the integrity of the monolayers by measuring trans-epithelial electrical resistance (TEER), and on differential gene expression relative to mock treated control using whole transcriptome sequencing technology. Our data establish complementary non-prebiotic modes of action for these frequently used dietary fibers.

Cells and Experimental Design
The microbiota-independent effects of GOS and FOS oligosaccharides were assessed in-vitro using human Caco-2 cells cultured as a confluent monolayer and exposed to treatment for 24 h. Each treatment was undertaken as three technical replications of three independent biological replicates represented by different cell passages. The human colon adenocarcinoma Caco-2 cell line was obtained from the American Type Tissue Collection (ATCC® HTB-37™, passage 18; lot 62381028; LGC Standards, Teddington, UK,) and propagated according to established methods [28]. In brief, cells were cultured in Dulbecco's modified Eagle medium (DMEM) and seeded at a density of 10 5 cells onto 1.12 cm 2 Transwell® permeable support with a 0.4 μm tissue culture treated polyester membrane placed in a 12-well plate (3460; Costar, Kennebunk, ME, USA). Caco-2 cells were maintained in a humidified atmosphere of 95% air and 5% CO2 at 37 °C. Fresh culture medium containing DMEM (D6546, Sigma, Gillingham, UK), sterile Fetal Bovine Serum 10% v/v (F9665, Sigma, Gillingham, UK), L-glutamine 2 mM (G7513, Sigma, Gillingham, UK), MEM non-essential amino acid solution supplement 1x (M7145, Sigma, Gillingham, UK), penicillin 100 U/mL and streptomycin 100 µ g/mL (P0781, Sigma, Gillingham, UK) was replenished every 2 to 3 days. Control mono-and di-saccharides-containing "GOS mock" and "FOS mock" treatment media were prepared as indicated in table 2. GOS and FOS treatment media were prepared by dissolving 2% (v/v) of the oligosaccharides syrups in FBS-free and antibiotic-free cell culture DMEM. All treatment culture media were free from serum and antibiotics to eliminate any interference from extraneous molecules, proteins or hormones.
After 18 days of culture, a confluent monolayer was obtained with a TEER exceeding 400 Ωcm 2 as measured by an EVOM voltohmmeter (World Precision Instruments, Sarasota, FL, USA). Independently, replicates of cells from passages ranging 22 to 26 were gently washed 3 times (at hourly intervals) with warm nonsupplemented (serum-free, antibiotic-free) DMEM culture medium and finally replenished from apical and basolateral sides of the Transwell® inserts with oligosaccharides or mock treatment media for 24h exposure before being harvested for total RNA extraction.

Integrity of Tight Junction Dynamics
Development of tight junction formation was monitored at regular intervals during the monolayers expansion by measuring the TEER across each of the growth inserts with an EVOM voltohmmeter (World Precision Instruments, Sarasota, FL, USA) connected to a pair of chopstick electrodes. Ohmic resistance of a blank (culture insert with respective medium but without cells) was measured in parallel. To obtain the sample resistance, the blank value was subtracted from the total resistance of the sample-containing cells. Final unit area resistance (Ωcm 2 ) was calculated by multiplying the sample resistance by the effective area of the membrane (1.12 cm 2 for 12-wells Transwell® inserts). To assess the effect of oligosaccharides on the Caco-2 cell monolayers, the electrical resistance across each cell monolayer was measured prior to exposure (T0) and after 24 hours of exposure (T24h) to the treatment media.

RNA Extraction and Qualification
Cells were gently washed in situ twice with ice cold sterile PBS followed by addition of 500 μl of TRIzol™ reagent (ThermoFisher Scientific, Paisley, UK) ready for harvesting by directly scraping the cells from the culture surface. Harvested cells were transferred to a 2 mL Eppendorf® safe-lock microcentrifuge tube for total RNA extraction according to TRIzol™ reagent manufacturer's protocol. Purified RNAs were quality-assessed using Agilent TapeStation 2200 system (G2964AA; Agilent, Stockport, UK) with RNA ScreenTape Assay reagents (5067-5576; 5067-5577; 5067-5578 Agilent, Stockport, UK). The RNA preparations had RIN values within the range 8.2 to 9.7 (mean = 9.24), which were subsequently aliquoted to limit freezethaw cycles and stored in −80 °C until cDNA library preparation.

Library Preparation for mRNA Sequencing
RNA concentrations were measured using Qubit Fluorometer and Qubit RNA BR Assay Kit (Q10211; ThermoFisher Scientific, Paisley, UK). A Biomek 4000 Automated Workstation (Beckman Coulter, High Wycombe, UK) was used to prepare sequencing libraries. mRNA was purified from 1 µ g of total RNA using the NEBNext Poly(A) mRNA Magnetic Isolation Module (E7490; New England Biolabs, Ipswich, UK). Indexed sequencing libraries were then prepared using the NEBNext Ultra directional RNA Library Preparation Kit for Illumina (E7760; New England Biolabs, Ipswich, UK) and NEBNext Multiplex Oligos for Illumina, Index Primers Sets 1 and 2 (E7335 and E7500; New England Biolabs, Ipswich, UK). Libraries were quantified using Qubit Fluorometer and Qubit dsDNA HS Kit (Q32854; ThermoFisher Scientific, Paisley, UK). Library fragment-length distributions were analyzed using the Agilent TapeStation 4200 with the Agilent High Sensitivity D1000 ScreenTape Assay (5067-5584 and 5067-5585; Agilent, Stockport, UK). Libraries were pooled in equimolar amounts and final library quantification performed using KAPA Library Quantification Kit for Illumina (KK4824; Roche, Pleasanton, CA, USA). The library pool was sequenced on an Illlumina NextSeq500 over two High Output 150 cycle kits (FC-404-2002; Illumina, San Diego, CA, USA), to generate over 40 million pairs of 75-bp paired-end reads per sample.

Quality Control and Mapping Analyses
High-throughput sequence reads were quality checked with FastQC (version 0.11.3). Raw reads in FASTQ format were subsequently processed using the Trim Sequences Module in CLC Genomics Workbench 12.03 (Qiagen, Aarhus, Denmark). In this step, clean data were obtained by removing reads from the raw data that contained adapter sequences and those of low-quality, all downstream analyses were performed with high-quality clean data. Using the RNA-Seq Analysis Module in CLC Genomics Workbench 12.03, clean trimmed reads were paired and mapped to Homo sapiens Hg38 (GRCh38 reference genome available in ENSEMBL) thus generating normalized counts of gene and transcript hits as RPKM (reads per kilobase of exon model per million mapped reads) and TPM (transcripts per million). Normalized counts were treated to compute P-values and fold changes (FC) for downstream differential expression analysis. Using the RNA-Seq Analysis Module from CLC Genomics Workbench 12.03, differential expression analysis between oligosaccharides-treated group and mock-treated group was determined using a negative binomial distribution model with the resulting P-values being adjusted to p-adj values using the approach of Benjamini and Hochberg [29] for controlling the false discovery rate (FDR). Genes with absolute FC ≥ 1.5 and p-adj value < 0.05 were designated as significantly differentially expressed (DE). Volcano plots of genes that were differentially expressed were developed using CLC workbench 12.03. Principal Component Analysis (PCA) was carried out using the web-based tool ClustVis [30]. Heatmap diagrams were developed using the Morpheus web-based tool [31] where rows and columns were clustered using Pearson's correlation distance and average linkage. Venn diagram showing the overlap of differentially expressed genes between GOS and FOS treatments were calculated using the web-based tool InteractiVenn [32].

Functional Analysis of Differentially Expressed Genes
Functional analysis of significant differentially expressed genes was achieved by implementing bioinformatics tools on gene set enrichment and pathway databases analysis based on Gene Ontology (GO) [33,34] terms and the Kyoto Encyclopedia of Genes and Genomes (KEGG) [35] collection of databases. We investigated statistically significant gene sets modulated by GOS and FOS with | FC value | ≥ 1.5 and FDR p-adj. < 0 .05 by implementing hypergeometric distribution within NIPA.Enrichment.R script (v0.6.7.R, https://github.com/ADAC-UoN/NIPA) and applying a minimum of 2 genes enriched, a cut-off q-value = 0.05 for KEGG and GO bases and an enrichment analysis established on "hsa" (Homo sapiens) Ensembl database. Briefly, NIPA.Enrichment.R script executed the bioconductor biomaRt [36,37] package to integrate Ensembl sequence data with data analysis using the GO tool, the Generally Applicable Gene set Enrichment for pathway (GAGE) package performed gene set analysis [38] and Pathview completed pathway-based data integration and visualization [39].

Quantitative Real-Time PCR Validation
Relative gene expression of DE genes computed by RNA-seq was established by RT-qPCR according to the Minimum Information for publication of Quantitative real-time PCR Experiments (MIQE) guidelines [40]. Sequences of RT-qPCR primers (Eurofins Genomics, Ebersberg, Germany) were designed using National Center for Biotechnology Information (NCBI) primer designing tool (http://www.ncbi.nlm.nih.gov/tools/primer-blast/) [41] and are described in table 3. cDNA was synthesized using Invitrogen™ SuperScript™ II (18064014, ThermoFisher Scientific, Loughborough, UK) according to the random hexamers primer method previously described [42]. Quantitative RT-PCR was performed using the Applied Biosystems PowerUP™ SYBR™ Green PCR Master Mix (A25742, ThermoFisher Scientific, Loughborough, UK) on a Real-Time PCR LightCycler® 480 System (Roche, Pleasanton, CA, USA). Each sample was processed in triplicate and the expression of genes of interest was normalized to endogenous housekeeping genes GAPDH, PUM1 and ACTB. Amplification protocol comprised one denaturation cycle at 95 °C for 5 min (ramp rate 4.8 °C/s), forty amplification cycles at 95 °C for 15 s, 65 °C for 30 s, 72 °C for 30s (ramp rate 1.5 °C/s). Melting curves were assessed for each RT-qPCR reaction using LightCycler® 480 System software (v 1.5.1.62; Roche, Pleasanton, CA, USA). Resulting PCR product size was confirmed using Agilent TapeStation 2200 system (G2964AA; Agilent, Stockport, UK) with DNA D1000 ScreenTape Assay (5067-5582 and 5067-5583; Agilent, Stockport, UK). Primers sets were validated for specificity when melting curves analysis exhibited a single peak with Tm > 78 °C and amplified product size matched expected product size as established by the NCBI primer designing tool. Relative gene expression was calculated using the relative quantification method with amplification efficiency corrected calculation models [43,44]. Amplification rate and relative concentration was calculated using LightCycler® 480 System software (v 1.5.1.62; Roche, Pleasanton, CA, USA) based on a linear regression slope established by 2-fold dilution series of a pool of all RNA samples.

Statistical Analysis
TEER experimental results are expressed as mean ± standard deviation, the differences between groups were evaluated by one-way ANOVA using Genstat 19th edition (VSNI Ltd, Hemel Hempstead, UK), with p-value < 0.05 was considered statistically significant. Within CLC Genomics Workbench 12.03, differential gene expression statistical analysis was based on Baggerly's beta-binomial model [45] to account for between-library variability; two-sided p-values for the test are described as 'P-value'. P-values were subsequently corrected for false discovery rate FDR ('padj values') using the Benjamini and Hochberg's approach [29]; p-adj value < 0.05 was considered significant. For functional analysis, hypergeometric distribution model generated q-values [46]; a minimum of 2 genes enriched with q-values < 0.05 for KEGG and GO enrichment were used as the criteria for reporting.

Influence of GOS and FOS on Epithelial Monolayer Tight Junction Integrity
To assess the direct effect of dietary GOS and FOS, cells were exposed to tailored media containing the oligosaccharides. The presence of mono-and di-saccharides contained in GOS and FOS syrups used for treatment was accounted for in the control mock media for Nutrabiotic® GOS 64% by incorporating galactose, glucose and lactose and by inclusion of fructose, glucose and sucrose for Orafti® L95 FOS. The treatment of Caco-2 cell monolayers with GOS generated a significant increase in the TEER values ( Figure 1, TEER +33.62%, p = 0.00037) when compared to the mock-exposed cells. Similarly, FOS also induced greater TEER when compared the mock-exposed cells (Figure 1, TEER +28.68%, p = 0.054). The rise in trans-epithelial electrical resistance of the monolayer suggest an acceleration of the tight junction dynamics that indicates an improvement of the monolayer integrity under the influence of the oligosaccharide treatments. Polarized confluent Caco-2 monolayers were exposed for 24 h to GOS, FOS (2% v/v) or their respective mock media. GOS increased significantly the monolayer TEER (+33.62%, p < 0.001), while the TEER increase induced by FOS exposure was not statistically significant (+28.68%, p = 0.054). Data are expressed as means ± SD and were tested using ANOVA t-test, ***p < 0.001 (n = 3).

Transcriptome Sequencing and Functional Annotation
The effect of GOS and FOS preparations on the global transcriptome of Caco-2 cells monolayers was assessed by RNA-seq using the Illumina NextSeq sequencing platform. RNAs were harvested from three biological replicates for each treatment group and mRNA sequencing generated between 44 to 134 million paired-end reads of 75-bp per sample. An average of 88.67% good quality paired-end reads per sample were mapped to Homo sapiens Hg38 reference genome. The summary of RNA-seq data and mapping for each sample are presented in Table 4. Gene and transcript normalized hit counts were used to generate FC and FDR adjusted p-values to implement differential gene expression analysis. Table 4. Summary of sequencing reads mapped to the human reference genome. RNA sequencing was strand-specific and only reverse strand reads were mapped to the reference genome Homo sapiens Hg38. FOS and GOS are oligosaccharides treated groups, the control groups are referred as "mock"; each group is represented by 3 independent biological replicates (n = 3).

Analysis of Differential Gene Expression
Volcano plots were generated to visualize the distribution of differentially expressed genes following GOS and FOS exposure (Figures 2A and 3A respectively), which showed a limited number of genes modulated by FOS when compared to GOS. Principal component analysis of the normalized transcript counts (TPM) disclosed tight relationship between biological replicates within each treatment group such that oligosaccharides-exposed samples clustered distinctly from their mock-treated counterparts. As a result, the first two principal components taken together explained 91.9% and 87.8% of the variability among the samples exposed to GOS and FOS when compared to mock treated cells ( Figures 2B and 3B respectively). We further examined the distribution of differentially expressed genes exhibiting FDR p-adj < 0.05 and |FC| ≥ 1.5 using heatmaps, which are presented in Figures 2C and 3C for the GOS and FOS treatments respectively. Following exposure to GOS, 89 genes were differentially expressed according to the criteria, of which 53 were up-regulated and 36 down-regulated when compared to control mock-treated cells ( Figure 2C and Table S1). Whereas for FOS, the total number of genes fulfilling the criteria for differential expression was limited to 12, comprising of eight upregulated and four down-regulated genes ( Figure 3C and Table S1). GOS and FOS treated cells shared five differentially expressed genes as presented in the Venn diagram shown in Figure 3D: TMEFF1, GPX2, SLC5A3, SULT1A3, AL139011.2. Three of these genes were up-regulated upon exposure to GOS and FOS treatments (GPX2, SLC5A3 and AL139011.2) and two genes showed the opposite modulation such that TMEFF1 transcripts were reduced upon GOS exposure whilst increased with FOS treatment, and SULT1A3 transcription was increased with GOS but reduced with FOS exposure. and y show principal components 1 (PC1) and 2 (PC2) that explain 77.1% and 14.8% of the total variance respectively. Prediction ellipses indicate probability 0.95 that a new observation from the same group will fall inside the ellipse; N = 6 data points. C. Heatmap displaying differentially expressed genes for GOS and mock treatments. Rows and columns are clustered using Pearson's correlation distance and average linkage. Data are presented as normalized transcripts count (TPM) with |FC| ≥ 1.5 and q-value < 0.05.

Validation of Differentially Expressed Genes by RT-qPCR
RNA-seq data was independently validated by RT-qPCR measurement of the transcript levels of genes identified as differentially expressed (RT-qPCR raw data presented in Table S2). Results revealed strong correlation between the RNA-seq and RT-qPCR differential expression data for the GOS and FOS treatments (

Functional Analysis of Differentially Expressed Genes
Biological attribute identification and functional interpretation of our data was achieved by implementing bioinformatics tools on gene set enrichment and pathway database analyses based on Gene Ontology (GO) terms and the Kyoto Encyclopedia of Genes and Genomes (KEGG) collection of databases. The KEGG knowledge base links gene catalogues to functional information such as metabolism, membrane transport, signal transduction and cell cycle in pathway assemblies. GO resource provides computable knowledge for biological functions of genes and gene products structured within biological process, molecular function and cellular component classes.
GO classification of the 89 differentially expressed (DE) genes associated with GOS exposure placed the up-regulated genes within 10 GO terms ( Figure 5) describing metabolic processes (daunorubicin, doxorubicin and progesterone metabolic processes, alkaline phosphatase activity and digestion), cell membrane transport (transepithelial chloride transport, amino acid transmembrane transport) and tissue homeostasis. KEGG analysis of the gene sets (Table 5) indicated the upregulated gene-associated pathways were enriched for metabolic processes featuring folate biosynthesis (hsa00790), thiamine metabolism (hsa00730), thyroid hormone synthesis (hsa04918), protein digestion and absorption (hsa04974), membrane transport with salivary (hsa04970) and pancreatic (hsa04972) secretion. The GOS down-regulated genes were placed in 5 GO terms ( Figure 5) that identified processes related to sulfation, 3'-phosphoadenosine-5'-phosphosulfate sulfotransferase, glutathione derivative biosynthesis, glutathione peroxidase activity and triglyceride catabolism. The down-regulated GO terms were further identified in KEGG (Table 5) as reduced glutathione metabolism (hsa00480), chemical carcinogenesis (hsa05204), drug and xenobiotics metabolism -cytochrome P450 (hsa00982 and hsa00980). A further seven DE genes were enriched in a unique GO cell compartment described as the collagen-containing extracellular matrix (GO:0062023, p-adj.= 0.018) involving three down-regulated and four up-regulated genes ( Table 6).

Discussion
Non-digestible oligosaccharides are prebiotics able to selectively stimulate the components of gut microbiota beneficial to host health [47]. A consequence of nondigestion is that they will reach the colon intact to exert any non-prebiotic effects on the intestinal epithelium. In this study, we addressed the hypothesis that oligosaccharides could directly affect the transcriptome of epithelial cells representative of the intestinal mucosa and have exploited RNA-seq methodology to describe the transcriptome of confluent Caco-2 cell monolayers exposed to the galacto-oligosaccharides and fructo-oligosaccharides.
The Caco-2 cell line model was chosen for its ability to form tight junctions, microvilli, thus expressing enzymes and transporters characteristic of enterocytes [48,49]. A number of in vitro studies have evaluated the impact of various prebiotics on the cell barrier integrity when subject to pathogen-induced barrier disruptions [50,51], inflammation inducing cytokines and lipopolysaccharides [52] or toxininduced challenge with the mycotoxin deoxynivalenol [53]. Here exposure to GOS (0.9% w/v) increased trans-epithelial electrical resistance of the monolayer suggesting GOS can improve the integrity of the tights junctions. FOS (1.3% w/v) exposed cells similarly increased trans-epithelial electrical resistance despite marginally failing to meet significance. We sought to understand these effects and any wider transcriptional responses to galacto-oligosaccharides and fructooligosaccharides by establishing differential gene expression using RNA-seq. Differential expression was determined in prebiotic treated Caco-2 cell monolayers with respect to tailored mock exposure to the saccharides present in the GOS and FOS preparations. As the output from these studies we opted to use the threshold of FDR adjusted probability p-adj < 0.05 and an absolute fold change ≥ 1.5 to identify differential expressed genes as basis for in silico approaches to assess the functional consequences. A total of 89 DE genes were identified between GOS and mock treated cells, whereas for FOS exposure the number was reduced to 12 DE genes. The RNAseq data was validated by RT-qPCR using specific primer sets. Strong correlations were observed between the fold-change determined for the two methods with R 2 = 0.96 for GOS and R 2 = 0.92 for FOS treatments. The DE data were used for GO term enrichment and KEGG analyses to identify transcriptional pathways responding to prebiotic oligosaccharide exposure.
FOS treatment resulted in 12 DE genes comprising of four down-regulated and eight up-regulated. GO terms identified the up-regulated genes SULT1A3, CYP1A1 and PLA2G4B as related to flavonoid metabolism, IRE1-mediated unfolded protein response and steroid metabolism. Interrogation of the KEGG database confirmed these three genes are functional in ovarian steroidogenesis and arachidonic acid metabolism (FDR q-val < 0.001) and general metabolic pathways (FDR q-val = 0.013). CYP1A1 encodes a Cytochrome P450 monooxygenase exhibiting high catalytic activity for steroid hormones [54] and hydroxylation of certain polyunsaturated fatty acids (PUFA) [55], whereas PLA2G4B encodes phospholipase A2 that preferentially catalyzes hydrolysis of arachidonoyl phospholipids to release lysophospholipids and fatty acids, with a preference for arachidonoyl phospholipids. The products of these genes also form a nexus with the phosphatidylinositol 4-kinases products of AC104662.2 and PI4K2B around cell signaling, lipid and membrane trafficking pathways [56]. In the absence of direct evidence for any transcriptional changes of the genes encoding tight junction proteins to affect the increases in TEER observed, it is conceivable that FOS exposure results in compositional changes in membrane sterols and phospholipids that affect the physical properties of the cell membrane, altering rigidity/fluidity leading to stabilization or disruption [57]. This contrasts with the direct effect, although not attributed to oligosaccharides, reported on the expression of CLDN4, encoding the epithelial cell tight junction protein Claudin 4, in Caco-2 cells exposed to complex food homogenates of apple, broccoli and button mushroom [58].
GOS treatment led to 53 up-regulated genes that associated with 10 GO terms. Genes AKR1C1, AKR1C2 and AKR1B1 were related to the terms daunorubicin and doxorubicin metabolic process, cellular response to jasmonic acid stimulus, progesterone metabolic process and digestion (FDR q-val < 0.005). These genes encode enzymes harboring oxidoreductase and aldo-keto reductase activities that target steroids and prostaglandins with varying levels of bile acid-binding affinity. Association to daunorubicin (DNR) and doxorubicin (DOX) metabolic processes may be a consequence of database bias to studies of cancer cells. In absence of DNR and DOX treatments, the implication of the GO enrichment association to such pathways remains unclear. Progesterone metabolic process and cellular response to jasmonic acid stimulus relate also to upregulation of AKR1C1 and AKR1C2 genes due to their catalytic activity related to androgens inactivation (progesterone and testosterone respectively). Despite its plant origin, jasmonic acid is an oxylipin signaling molecule often considered a structural and functional analogue to prostaglandins in animals. Mammalian eicosanoids and oxylipins such as prostaglandins and leukotrienes are active signaling molecules derived from oxygenated PUFAs and are potent regulators of host immune responses. Notably, mammalian prostaglandins and leukotrienes (products of arachidonic acid metabolism) can polarize macrophages, modulate T helper cell immune responses, stimulate chemokine production, phagocytosis, lymphocyte proliferation and leukocyte chemotaxis [59].
GO also associated AKR1C1 and AKR1C2 with CAPN8 within digestion processes. CAPN8 encodes tissue-specific calpain 8a (cytosolic calcium-dependent cysteine protease), which in mice is reported to be expressed specifically in the mucus-secreting pit cells of the gastric mucosa as well as in a subset of goblet cells in the small intestine [60]. Calpains are also implicated in membrane trafficking processes due to their localization in endoplasmic reticulum, Golgi and lipid rafts [61,62]. The KEGG protein digestion and absorption pathway (hsa04974) identified the genes SLC7A8, COL12A1 and FXYD2, which respectively encode L-type amino acid transporter 2, collagen type XII alpha 1 chain and Na/K-transporting ATPase subunit gamma suggesting specific amino acid transport and focused structural protein remodeling. The KEGG pancreatic secretion (hsa04972) and salivary secretion (hsa04970) pathways featured the up-regulated genes SLC12A2 encoding a basolateral Na-K-Cl symporter and ADCY1 encoding Ca 2+ /calmodulin-activated adenylyl cyclase. These activities are linked as Ca 2+ and cAMP can regulate intestinal Na-K-Cl symporter activity [63]. Related to this the GO terms transepithelial chloride transport (FDR q-val < 0.001) and amino acid transport (FDR q-val < 0.05) identified genes for ion channel and solute carrier proteins (BEST1, SLC12A2, SLC6A20, SLC7A8 and SLC6A12). Several solute carrier encoding gene were up regulated: SLC7A8, SLC6A20, SLC6A12, SLCO4A1, SLC12A2, SLC7A8, SLC6A20, SLC5A3, SLC6A12 and BEST1; that feature in Reactome pathway analyses related to amino acid transport across the plasma membrane and the transport of small molecules. Two bioprocesses relevant to the digestion pathway were enriched in KEGG folate biosynthesis (hsa00790) and thiamine metabolism (hsa00730) pathways. Genes ALPP and ALPG encoding intestinal and placental-like alkaline phosphatase respectively were identified with AKR1B1 encoding aldo-keto reductase. AKR1B1 in association with SOX9 was identified as enriched constituents of the GO tissue homeostasis pathway. SOX9 encodes SRY-Box Transcription Factor 9 that binds to the COL2A1 promoter to activate cartilaginous tissue-specific Collagen Type II Alpha 1 Chain protein expression, which has an essential role in skeletal development and contributes to the ability of cartilage to resist compressive forces. Collectively these data are consistent with the hypothesis that GOS elicits energy-dependent transmembrane trafficking of solutes in Caco-2 cell monolayers accompanied by collagen and cytoplasmic membrane remodeling that may be linked to the observed increases in TEER.
Of 36 DE genes found down-regulated following GOS exposure, seven DE genes were annotated in enriched GO bioprocess terms associated to sulfation & 3'phosphoadenosine 5'-phosphosulfate metabolic process (genes SULT1E1 and SULT1C2), glutathione peroxidase activity and derivative biosynthetic process (genes GSTA2, GSTA1 and ALOX5AP), and triglyceride catabolic process (genes APOC3 and FABP2). All these genes encode proteins related to fatty acid and arachidonic acid metabolism. Genes SULT1E1 and SULT1C2, encode sulfotransferases responsible for sulphate conjugation of many hormones (such as estrogens), neurotransmitters and xenobiotic compounds. GSTA2 and GSTA1 encode alpha-glutathione S-transferases (GSTs) involved in the synthesis of prostaglandins (PGs) and leukotrienes thus responsible for the metabolization of electrophilic compounds including carcinogens, therapeutic drugs and the byproducts of oxidative stress. Notably, through eicosanoid metabolites GSTs can play a central role as regulators of transcription factor peroxisome proliferator-activated receptor γ (PPARγ) [64,65]. In mammals, overexpression of GSTs in tumor cells has been implicated with resistance to various anti-cancer agents, chemical carcinogens and in microbes, plants and mammals expression of GSTs are up-regulated by exposure to pro-oxidants, thus suggesting that induction of GST is an evolutionary conserved response of cells to oxidative stress [66]. ALOX5AP encodes arachidonate 5-lipoxygenase activating protein required for leukotriene biosynthesis by 5lipoxygenase (ALOX5). The products of LOX metabolism either represent biologically active eicosanoid metabolites such as hydroperoxy-eicosatetraenoic acids (HPETEs) or give rise to the production of reactive oxygen species (ROS). A functional role has also been attributed to lipoxygenase-catalyzed arachidonic and linoleic acid metabolism in cancer development [67,68]. APOC3 encodes apolipoprotein C3, a protein component of the triglyceride (TG)-rich lipoproteins (TRLs) and FABP2 encodes intestinal-type Fatty Acid-Binding Protein -an intracellular fatty acid-binding protein that participates in the uptake, metabolism and transport of long-chain fatty acids. This lends further support to the reconfiguration of fatty acid metabolism in response to GOS treatment whilst maintaining homeostasis.
GO cell compartment analysis identified the collagen-containing extracellular matrix based on the downregulation of APOC3, ZG16 and ADAMTS9 and upregulation of COL12A1, S100A4, S100A6 and SERPINA5. Gene APOC3 encodes a liver and small intestine-secreted apolipoprotein C-III (apoC-III) found with triglyceride-rich lipoproteins, such that increases in apoC-III levels induces the development of hypertriglyceridemia and overexpression can contribute to coronary heart disease in humans [69,70]. Following exposure to GOS, APOC3 was down-regulated thus suggesting that a reduction in apoC-III levels could be beneficial for triglyceridemia. However, ZG16 and ADAMTS9, which have been ascribed roles in the control of cancer development, were also downregulated. Zymogen granule protein 16 (ZG16) is one of the most significantly downregulated genes in colorectal cancer (CRC) and harbors a jacalin-like lectin domain. Lectins are carbohydrate-binding proteins able to detect subtle differences between complex carbohydrate structures thus recognizing specific sugars to carry out various functions including cell attachment, migration and invasion. Amongst these, some galectins share an affinity for simple β-galactoside moieties [71]. ZG16 may be an important component of the protective mucus layer, which helps separate host epithelium from commensal bacteria in the colon [72]. The high similarity of ZG16 to jacalin suggests that the human homologue may play an important role in colon cancer immunity. The loss of ZG16 associated with CRC development has led to the hypothesis that ZG16 reduction may disrupt well-organized bacterial surveillance systems to facilitate bacterial invasion of host tissues and cause local inflammatory changes, which may constitute an increased risk for the development of cancer [73]. ADAMTS9 encodes a member of the ADAMTS (a disintegrin and metalloproteinase with thrombospondin motifs) protein family, a secreted mammalian metalloprotease that localizes to the cell-surface and/or extracellular matrix. The extracellular matrix is continuously remodeled by coordinated biosynthesis and proteolysis of its components. Active proteolytic ADAMTS9 has been shown to interact with fibronectin and disrupt fibril networks [74]. Furthermore, it has been suggested ADAMTS9 contributes to the suppression of tumorigenesis by decreasing cell proliferation, inducing cell apoptosis and inhibiting angiogenesis through regulating AKT/mTOR signaling pathway, whilst methylation of ADAMTS9 genes is associated with poor survival of gastric cancer patients [75]. These data must be considered in the context that Caco-2 cells represent colonic cancer cells, and that the regulatory differences favoring cell proliferation could be a consequence of this.
Collagen-containing extracellular matrix enriched genes S100A4 and S100A6 were found to be upregulated following treatment with GOS. Proteins S100A4 and S100A6 are members of the S100 calcium binding protein family that exert biological function via interactions with target proteins. Many kinds of cells including fibroblasts, immune cells, and cancer cells can produce S100A4, which are released into the extracellular space in response to various stimuli such as activated normal T cell expression and secreted factors (RANTES) produced by the tumor cells. S100A4 could be involved in the regulation of cell proliferation and differentiation, apoptosis, Ca 2+ homeostasis, and energy metabolism [76]. Sun et al. demonstrated mice deficient in S100A4 exhibited impaired humoral and cellular immunity after mucosal immunization using Cholera toxin as adjuvant, revealing a reduced production of Th1, Th2, and Th17 cytokines [77]. Protein S100A6 has been described in a limited number of cell types in adult normal tissues and in several tumor cell types. As an intracellular protein, S100A6 has been implicated in the regulation of several cellular functions, such as proliferation, apoptosis, cytoskeleton dynamics, and cellular response to different stress factors. S100A6 can be secreted/released by certain cell types which points to extracellular effects of the protein such as antimicrobial activity [78]. Similar upregulation was observed for gene SERPINA5 encoding Protein C inhibitor (PCI), a serine protease inhibitor of serpin type that is found in most tissues and fluids, including blood plasma, seminal plasma and the urine of humans [79]. As an antimicrobial agent, PCI has the ability to disrupt the bacterial cell wall to causing death by interacting with lipid membranes leading to permeabilization of bacterial pathogen [80]. PCI also inhibits proteases of the blood coagulation and fibrinolysis system whilst in cancer cells it suppresses tumor invasion by inhibiting urokinase-type plasminogen activator and inhibits tumor growth and metastasis, which are independent of its protease-inhibitory activity [81]. Lastly, COL12A1 encoded collagen XII is a member of the family of fibrilassociated collagens with an interrupted triple helix (FACIT) structure. Mutations in the collagen XII gene associated with extracellular matrix-related myopathy, as the protein functions to preserve muscle and bone architecture through collagen fibril organization [82,83].

Conclusion
In this study, the direct effects of galacto-oligosaccharides and fructooligosaccharides on colonic epithelial cells were assessed through changes in monolayer transepithelial resistance and transcriptome analysis. Results suggest galacto-oligosaccharides have the potential to directly increase integrity of the epithelial barrier. Transcriptome data suggests TEER increased independent of specific tight junction protein expression but could be due to sterol/fatty acid compositional changes associated with transporter dependent remodeling of the cell membrane. Specific candidate pathways group the genes involved in digestion and transepithelial transport, which contribute to intestinal cell integrity and function.
Carbohydrates have enormous potential to encode biological information. It is conceivable that galacto-oligosaccharides and fructo-oligosaccharides can crosstalk with cells using lectin molecules as the interface. Typically, lectins and their complimentary carbohydrate are located on the surfaces of opposing cells, which can be of the same type or different types. Such interactions are required for cell differentiation, development and pathological states. These results highlight concerted effects on transmembrane trafficking, differences in xenobiotic biotransformation, and the production of antimicrobial agents.
GOS and FOS treatments shared relatively few differentially expressed genes suggesting they have different modes of action. Our strategy has produced a comprehensive database of gene expression profiles of Caco-2 cell monolayers exposed to oligosaccharide food ingredients, allowing further work to link gene expression signatures of cultured cells to their mode of action, and thus potentially facilitating product choice in human or animal intervention studies.

Supplementary Materials:
The following are available online at www.mdpi.com/xxx/s1, Supplementary Table S1: GOS and FOS gene expression browser, Table S2: RT-qPCR raw data.
Author Contributions: Conceptualization, NMF and IFC; methodology, GFL and IFC; software, GFL; validation, GFL and IFC; formal analysis, GFL and IFC; investigation, GFL; resources, NMF and IFC; data curation, GFL and IFC; writing-original draft preparation, GFL; writing-review and editing, GFL, NMF and IFC; visualization, GFL; supervision, IFC; project administration, GFL, NMF and IFC; funding acquisition, NMF and IFC. All authors have read and agreed to the published version of the manuscript.
Funding: This research was funded by Saputo Dairy UK.