Differential Gene Expression Profiling of Orbital Adipose Tissue in Thyroid Orbitopathy

Abstract

Purpose: We aimed to determine differentially expressed genes relevant to orbital inflammation and orbital fat expansion in thyroid orbitopathy (TO) using microarray gene profiling in a case-control study.

Methods: Human orbital adipose samples were obtained from individuals with active TO (n = 12), inactive TO (n = 21), and normal controls (n = 21). Gene expression profiles were examined using microarray analysis and were compared between active and inactive TO, and between active TO and normal controls. Top ranked differentially expressed genes were validated by real-time RT-PCR in an additional eight active TO, 13 inactive TO, and 11 normal controls and correlated with gene set enrichment analysis (GSEA) and molecular pathways analysis.

Results: Seven hundred twenty-one probes (683 genes) and 806 probes (735 genes) were significantly differentially expressed in comparing active to inactive TO and in comparing active TO to healthy controls, respectively. All selected genes were confirmed to be differentially expressed by real-time RT-PCR. Multiple top ranked genes in active versus inactive TO comparison are overrepresented by immune and inflammatory response genes. They include defensins (DEFA1, DEFA1B, DEFA3), which were overexpressed by 3.05- to 4.14-fold and TIMD4 by 4.20-fold. Markers for adipogenesis were overexpressed including SCD, FADS1, and SCDP1. Gene set enrichment analysis revealed dysregulation of epigenetic signatures, T-cell activation, Th1 differentiation, defensin pathway, cell adhesion, cytoskeleton organization, apoptosis, cell cycling, and lipid metabolism in active TO.

Conclusions: Active TO is characterized by upregulation of genes involved in cell-mediated immune, innate immune, and inflammatory response and enhanced orbital adipogenesis. TIMD4, DEFA1, DEFA1B, and DEFA3 genes may be involved in the innate immune-mediated orbital inflammation in TO. Epigenetic mechanisms may play a role in the pathogenesis of TO.

Thyroid orbitopathy (TO) is an autoimmune disorder of the orbit related to Graves' disease. Thyroid orbitopathy occurs in 25% to 50% of patients with Graves' disease, and 3% to 5% present with potentially blinding complications such as optic nerve compression or exposure keratopathy.¹ Extraocular muscle enlargement from glycosaminoglycan deposition and interstitial swelling, connective tissue inflammation, and orbital fat expansion are part of the observed pathologic process, which leads to periocular edema, proptosis, epiphora, retro-orbital pain, and diplopia.^2–5 The molecular mechanisms involved in orbital fat tissue expansion remains incompletely understood. Previous microarray studies suggested peroxisome proliferator activated receptor gamma, secreted frizzled-related proteins, adipocyte-related immediate early genes, and lysosome-related genes may play a role in orbital adipogenesis.^6–9 

Expression microarray profiling uses high throughput technology to generate extensive gene expression data from small amounts of tissue. It is a rapid and highly effective way to discover differentially expressed genes with high sensitivity, thus advancing the understanding of pathogenesis without biases from prior knowledge of the disease. We aimed to determine differentially expressed genes that may be involved in stimulating orbital inflammation and orbital adipose tissue expansion in TO using microarray gene expression profiling in a case-control study to improve understanding of the molecular processes involved in this disease. 

Methods

Tissue Samples

Human orbital fat samples were obtained during orbital decompression (extraconal orbital fat) and upper lid surgeries (preaponeurotic orbital fat) in patients with TO. Patients with TO were subclassified as active (n = 12) or inactive (n = 21). Normal controls (n = 21) were patients without autoimmune thyroid disease with adipose tissue harvested from corresponding anatomical locations at the time of unrelated orbital and lid surgeries. Thyroid orbitopathy status was determined by an ophthalmologist according to vision, inflammation, strabismus, appearance (VISA) classification.¹⁰ Active TO is defined by the presence of signs of orbital inflammation; inactive TO was defined by the absence of these inflammatory signs. Activity of thyroid eye disease was defined by inflammatory index score; the maximum inflammatory index score is 8, where chemosis could score up to 2, conjunctival redness up to 1, lid erythema up to 1, lid edema up to 2, and retrobulbar pain up to 2 points. Seventy-five percent of active TO and 76% of inactive TO are euthyroid, 25% of TO patients had subclinical hyper- or hypothyroidism while on medical treatment, but there was no statistically significant difference in the means of TSH and T4 levels. The research was conducted in compliance with the Declaration of Helsinki and was approved by human research ethics committee of the Royal Victorian Eye and Ear Hospital, Melbourne, Australia. 

Following excision, orbital adipose tissue was immediately immersed in RNAlater solution (Ambion, Austin, TX, USA). The specimens were stored at 4°C for at least 24 hours before long-term storage at −80°C. 

RNA Extraction and Hybridization to Microarrays

Total RNA was extracted using RNeasy microarray tissue Mini Kit 50 (Qiagen, Valencia, CA, USA) according to the manufacturer's protocol. The quality and quantity of RNA was ascertained on a Bioanalyzer 2100 (Agilent Technologies, Inc., Santa Clara, CA, USA) using NanoChip according to manufacturer's protocol. The mean RNA Integrity Number (RIN) for the RNA samples was 6.92 ± 0.7 SD. RNA samples were processed and hybridized to microarrays in two batches 1 year apart at the Australian Genome Research Facility using Illumina TotalPrep RNA Amplification Kit (Life Technologies, Carlsbad, CA, USA) protocol. A total of 750 ng amplified RNA was prepared for hybridization to HumanHT-12 v4 Expression Beadchip (Illumina, San Diego, CA, USA) by preparing a probe cocktail (biotin labeled cRNA at 0.05 μg/μL) that includes GEX-HYB Hybridization Buffer supplied with the beadchip. A total hybridization volume of 15 μL for each sample was loaded into a single array on the beadchip, then hybridized at 58°C for 16 hours and coupled with Cy3 for scanning in the Illumina iScan Reader. Raw intensity values for both regular and control probes were exported to text files using Illumina GenomeStudio software. 

Microarray Data Analysis

Bioinformatics analysis was conducted using the limma software package.¹¹ The raw intensity values were background corrected and normalized using the neqc function, which performs “normexp” background correction and quantile normalization using parameters estimated from the control probes.¹² One probe (targeting long intergenic nonprotein coding RNA 1239) was removed because it was available only in the batch 2 files. The proportion of probes expressed in each sample was estimated using the propexpr function, which compares regular probes to negative controls.¹³ The average proportion of expressed probes was found to be 53.6%, so the top 53.6% probes with highest average log normalized expression values were retained for subsequent analysis. This left 25,796 probes. The Illumina manifest file (HumanHT-12_V4_0_R2_15002873_B) was used to associate an Entrez Gene identifier with each probe. Other gene annotation was obtained from the National Center for Biotechnology Information Homo sapiens gene information file downloaded on November 21, 2013. 

Linear models were used to test for expression differences between disease categories while adjusting for sex differences and for differences between the hybridization batches. The correlation between repeated measurements from the same patient was estimated as 0.4 and incorporated into the linear models.¹⁴ Empirical array quality weights were estimated to allow for differences in quality between the RNA samples.¹⁵ Differential expression between active TO, inactive TO, and normal tissue was assessed using empirical Bayes moderated t-statistics,¹⁶ allowing for an intensity-dependent trend in the standard errors. The false discovery rate (FDR) was controlled at less than 0.05 using the Benjamini and Hochberg method. Genes associated with the inflammatory index score were identified the same way, except that the quantitative index score was used in the linear model instead of the disease categories. 

Molecular signature enrichment analyses were conducted using limma's mroast function,¹⁷ which uses residual rotation to conduct gene set tests. Roast was run with 99,999 rotations, using the same linear model settings as for the differential expression analysis including batch correction, patient correlations, and quality weights. Gene sets representing expression signatures were downloaded from the c2-curated collection of Version 4.0 of the Molecular Signatures Database.¹⁸ 

The lists of differentially expressed genes were uploaded to DAVID 6.7 webtool,^19,20 to test for overrepresentation of pathways from the Biological Biochemical Image Database (BBID), BioCarta, and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases. 

Quantitative Reverse Transcription–Polymerase Chain Reaction (qRT-PCR) Validation

A total of eight active, 13 inactive TO, and 11 normal control orbital adipose samples, with RNA extracted as above, were used to confirm findings. This included three active TO, five inactive TO, and three normal controls also used in the microarray study. All RNA samples underwent in-solution DNase digestion using RNAse-Free DNase Set (Qiagen) followed by on-column RNA cleanup using RNeasy Mini Kit (Qiagen). The mean RIN for these RNA samples was 6.72 ± 0.35 SD. 

Quantitative RT-PCR was performed in a 384-well plate with 19 gene targets chosen from the top ranked and some lower ranked genes of interest. Expression levels of SREBF1, TIMD4, DEFA3, DEFA1B, DEFA1, FADS1, SCD, ELOVL6, DYNCII1, PLOD2, CHDH, LTB, CD3D, CAMP, EOMES, CD8A, CCL5, GZMA, and SLAMF6 were compared in duplicate. cDNA synthesis was performed using RT² First strand Kit, and cDNA amplification was performed using RT² SYBR Green qPCR Mastermix in an Applied Biosystems viiA 7 (Life Technologies, Carlsbad, CA, USA) thermal cycler according to RT² Profiler PCR Array protocol (Qiagen) according to MIQE guidelines.²¹ Fold change and P values were analyzed using the web-based resource GeneGlobe Data Analysis developed by SABioscience (http://www.sabiosciences.com/pcrarraydataanalysis.php; provided in the public domain by SABioscience, a Qiagen company). 

Results

The study was designed to determine the differences in genes expression profiles in TO by comparing diseased orbital adipose tissue when TO is in the active and inactive phases, and by comparing active TO with normal healthy tissues in a case-control study design. Orbital adipose tissue samples were harvested during orbital decompression and preaponeurotic orbital fat during eyelid surgeries in TO cases. Orbital fat from corresponding anatomical positions was harvested from unaffected healthy controls. TO status was determined by an ophthalmologist based on grading of vision, signs of orbital inflammation, strabismus and ocular motility, and severity of appearance according to VISA classification.¹⁰ Active TO is defined by presence of signs of orbital inflammation including periocular edema, lid erythema, conjunctival chemosis, and injection; inactive TO defined by the absence of these inflammatory signs. 

Differentially Expressed Genes Between Active and Inactive TO Identified by Microarray Analysis

The Illumina HumanHT-12 v4 contains 47,323 probes. Overall, the average proportion of expressed probes was 53.6%, hence the top 53.6% probes equivalent to 25,796 probes were retained for analysis. Seven hundred twenty-one probes (355 up and 366 downregulated) representing 683 annotated genes were significantly differentially expressed in the active TO versus inactive TO samples (FDR < 0.05). The top differentially expressed genes were dominated by overexpressed cell-mediated and innate immune response genes and adipogenesis markers (Table 1). Twenty of the 40 top differentially expressed probes represented immune and inflammatory response genes. Of these, 19 were upregulated and one was downregulated. TIMD4 was overexpressed by 4.20-fold and DEFA3 by 4.14-fold. DEFA1B was represented by three different probes and was upregulated by 3.05- to 4.04-fold and DEFA 1 by 3.95-fold. Genes overexpressed by 1.5- to 2-fold were CST7, CD247, NKG7, PTPRCAP, EOMES, GZMA, CYTL1, SLAMF6, KLRG1, CCL5, and CAMP. MASP1, a C4/C2 activating component of Ra-reactive factor, was significantly downregulated. Markers of adipogenesis were overexpressed by 2.36- to 6.12-fold including FADS1, SCD, and SCDP1. Lipid synthesis regulatory genes SREBF1 and DBI were upregulated by 1.40- and 1.52-fold, respectively. FZD9 was downregulated by 1.50-fold, whereas GPX3 was upregulated by 1.83-fold. 

Table 1

 

Top 40 Differentially Expressed Gene Probes in Orbital Adipose Tissue in Active Versus Inactive TO

 

Differentially Expressed Genes Between Active TO and Normal Controls Identified by Microarray Analysis

Of the 25,796 expressed probes, 806 probes (429 up- and 377 downregulated) representing 735 genes were differentially expressed in active TO compared to normal controls (FDR < 0.05). Ten of the 40 top ranked probes were from genes involved in adipogenesis, specifically unsaturated fatty acid and cholesterol synthesis including SREBF1, SCD, SCDP1, FADS1, PTPLB, and ACSS2; these genes were upregulated by 1.57- to 5.37-fold (Table 2). Differential expression of other genes involved in metabolism included upregulation of PDXK and downregulation of CHDH. Eight of the top 40 ranked differentially expressed probes were from genes involved in cytoskeleton organization, collagen binding, cellular adhesion, and migration. Upregulated genes were PLOD2, HOOK2, and SYNC with a fold difference of 1.45 to 1.91. Downregulated genes were FNDC5, DYNC1I1, TMSB15B, COLGALT2, and RADIL with fold change ranging from 1.51 to 2.10. Five immune and inflammatory response genes, TIMD4, KLRG1, DEFA3, PTGER3, CST7, were upregulated by 1.51 to 3.70-fold. Two genes coding protein with zinc binding domain (PDZRN4) and ubiquitin-protein ligase activity (TRIM9) were downregulated, and the specific functions of these genes have not been determined. 

Table 2

 

Top 40 Differentially Expressed Gene Probes in Orbital Adipose Tissue in Active TO Versus Normal Controls

 

Correlation of Differential Expressed Genes With Inflammatory Index Score

To explore the effects of more graduated changes in orbital inflammation, we directly correlated the expression of each gene with the individual inflammation index scores instead of categorizing patients into normal, active TO, and inactive TO. This analysis detected 188 probes positively correlated with inflammation and 162 negatively correlated (FDR < 0.05; Supplementary Table S1). Genes previously found to be upregulated in active TO were generally found to be positively correlated with inflammation, while those previously found to be downregulated in active TO were negatively correlated. Supplementary Table S2 shows correlation results for selected top ranking genes from the disease comparisons including all 19 genes selected for PCR validation. These results confirm negative correlations for CHDH and DYNC1I1 and positive correlations for the other genes. 

Validation of Differential Expression of Genes in Orbital Adipose Tissue Specimens

Seventeen and eight gene targets were chosen for validation from the active versus inactive and active versus normal comparisons, respectively. Quantitative RT-PCR was used to analyze the expression of these genes in eight active TO, 13 inactive TO, and 11 normal control orbital adipose samples. The qRT-PCR results were closely correlated with the microarray results, with all log 2-fold changes of comparable size and in the same direction as the microarray results for the same genes. This was true both for the active versus inactive TO comparison (Table 3) and for the active TO versus normal comparison (Table 4). Most chosen genes were confirmed to be significantly differentially expressed in the validation samples; even those genes that failed to achieve statistical significance still showed expression changes of similar size and direction to the microarray results. 

Table 3

 

Validation of Differential Express of Genes in Active Versus Inactive TO by qRT-PCR

 

Table 4

 

Validation of Differentially Expressed Genes in Active TO Versus Normal Healthy Control by qRT-PCR

 

Gene Set Enrichment Analysis (GSEA) of Molecular Signatures

To understand the active TO expression profile more deeply, we conducted a GSEA using curated signatures from the Molecular Signatures Database.¹⁸ Each signature was tested using rotation gene set tests (ROAST),¹⁷ which are able to adjust for sex differences, patient effects, and batch effects. 

A number of epigenetic signatures were found significantly dysregulated in active TO, both in comparing to inactive TO and normal control. Upregulated hypomethylation of gene clusters in immunodeficiency signaling, cytotoxic T-cell mediated apoptosis and T-cell receptor, and hypomethylation of gene cluster representing epigenetic signature of PML-RARA leukemia (Figueroa AML methylation cluster 5 and cluster 6) were noted in active compared to inactive TO. Conversely, active TO compared to normal showed negative correlation with unmethylated histone H3 in genes with high-CpG-density promoters (HCP) in embryonic fibroblast, downregulation of histone H3 dimethylation at K4 and trimethylation at K27 in brain, and downregulation of genes methylated aberrantly in colon cancer cells (Supplementary Tables S3, S4). 

Active TO versus inactive TO was characterized by upregulation of cell adhesion and cytoskeleton organization pathways including stathmin, cell-to-cell, E-cadherin nascent adhesion junction, and e-cadherin stabilization pathway. Molecular processes for activation of T cell (CSK, TOB1, T-cell receptor activation pathways), cytotoxic T cell (downstream CD8 T cell, cytotoxic T-cell pathway), T helper cell function, IL-12 dependent Th1 development, Th1 specific genes and natural killer cell response, and immunoregulation in lymphoid cells were upregulated in active versus inactive TO. Interestingly, CD40 signaling was negatively correlated with active TO. Gene sets controlling cell proliferation and cell cycling were downregulated as exemplified by downregulation of ERBB network, nuclear ERBB4 signaling, but genes involved in mitogen-activated protein kinases (MAPK) signaling, neuronal proliferation, and differentiation were upregulated. Apoptosis was regulated by upregulation of caspase cascade and upregulation of genes negatively correlated with telomere shortening and telomerase reverse transcriptase, respectively. Lymphogenesis associated genes were upregulated. Multiple oncogenic genes, tumor suppressor genes, and homeobox genes were dysregulated in active TO in the context of comparing with cancer gene expression profiles (Supplementary Table S3). 

Analysis of the active TO versus normal control expression changes revealed correlation with a number of signatures regulating lipid metabolism. These include genes involved in lipid digestion, mobilization and transport, integration of energy and metabolism, SREBF and peroxisome proliferator-activated receptor (PPAR) targets, the WNT noncanonical pathway, and biosynthesis of unsaturated fatty acids. Cytoskeleton control stathmin pathway was upregulated. Apoptosis pathways were upregulated with increased proportion of genes in caspase cascade, upregulation of proapoptotic genes regulated by telomerase reverse transcriptase, camptothecin, and methotrexate. Cell cycling was downregulated with downregulation of genes involved in activation of prereplicative complex and Rad3-related kinase (ATR), G2 checkpoint in mitosis, and mini-chromosome maintenance pathway, whereas cell differentiation was upregulated (osteoblast differentiation in response to phenylamil). Genes related to hypoxia were upregulated. Generally, cancer-related gene profiles are downregulated in active TO compared to normal (Supplementary Table S4). 

Pathway Analysis Using DAVID Functional Annotation Analysis

Functional annotation was examined more systematically for the entire set of differentially expressed genes by using the DAVID tool.^19,20 Analysis of differentially expressed genes for the active versus inactive TO comparison revealed overrepresentation of the following pathways: KEGG's regulation of actin cytoskeleton (17 genes, P = 2.9 × 10⁻³), adherens junction (eight genes, P = 0.017), glycolysis (seven genes, P = 0.017), MAPK signaling (17 genes, P = 0.022), glutathione metabolism (six genes, P = 0.028), biosynthesis of unsaturated fatty acids (four genes, P = 0.039), focal adhesion (13 genes, P = 0.045); BioCarta's T cytotoxic cell surface molecules (five genes, P = 1.7 × 10⁻³), Lck and Fyn tyrosine kinases in initiation of TCR activation (four genes, P = 0.013), T-helper cell surface molecules (four genes, P = 0.016), CTL-mediated immune response against target cells (four genes, P = 0.031), and IL-17 signaling pathway (four genes, P = 0.031). 

Functional annotation of differentially expressed genes in active TO versus normal controls highlighted overrepresented molecular pathways involved in KEGG's PPAR signaling (13 genes, P = 5.6 × 10⁻⁶), citrate cycle (seven genes, P = 7.6 × 10⁻⁴), biosynthesis of unsaturated fatty acids (six genes, P = 9.8 × 10⁻⁴), insulin signaling (13 genes, P = 3.6 × 10⁻³), focal adhesion (16 genes, P = 6.4 × 10⁻³), butanoate metabolism (six genes, P = 7.3 × 10⁻³), pentose phosphate (five genes, P = 0.012), glycerolipid metabolism (six genes, P = 0.023), and glutathione metabolism (six genes, P = 0.035). 

Pathways for regulation of actin cytoskeleton, adherens junction, MAPK signaling, Lck and Fyn tyrosine kinases in initiation of TCR activation, PPAR signaling, and insulin signaling pathways are illustrated in Supplementary Figures S1 through S6 with arrows pointing up for upregulated and down for downregulated gene probes. The corresponding Supplementary Tables S5 through S10 correlates dysregulated gene probes with the coded proteins shown on the illustrated pathways, and showed extent of gene probe dysregulation. 

Discussion

This differential gene expression profiling study utilizes the largest number to date of orbital adipose tissues in TO patients and normal subjects. This is the first gene profiling study to compare TO in early and late phases and uses high density Illumina microarray platform to derive meaningful findings. The high throughput technology combined with selection of phenotypically well-characterized TO patients generated extensive gene expression data, which allowed examination of differentially expressed genes with high sensitivity. The results were validated by qRT-PCR. Functional annotation of genes showed good correlation between different bioinformatics analysis. 

The principal finding of this study is that cell-mediated and innate immune response genes and adipogenesis markers are concurrently overexpressed in active TO, and epigenetic pathways may be involved in the pathogenesis of TO. Among these, TIMD4, DEFA3, DEFA 1, and DEFA1B were highly overexpressed. T-cell immunoglobulin and mucin domain containing 4 (TIMD4, also known as TIM4) is a cell membrane receptor exclusively on antigen-presenting cells including dendritic cells and macrophages that recognizes phosphatidyl serine, a specific marker of apoptosis, and via its immunoglobulin domain has a role in phagocytosis of apoptotic cells and in immune tolerance.^22,23 Blockade of TIMD4 results in increased numbers of antigen specific T cells following the peak of immune response. Conversely, overexpression of TIMD4 on antigen-presenting cells resulted in decreased antigen-specific T cells and decreased secondary T-cell response in TIMD4 transgenic mice.²⁴ Overexpression of TIMD4 in peripheral blood mononuclear cells was noted in patients with systemic lupus erythematosus (SLE), especially during the active phase of the disease compared with healthy controls. Furthermore, TIMD4 mRNA expression positively correlated with TIM1 and TNF-α in SLE but not in healthy controls.²⁵ 

TIM1 expression in polarized human Th2 cells is known to be upregulated in allergic rhinitis patients compared with healthy controls. TIM4 on dendritic cells played a critical role in maintaining skewed Th2 response by interacting with TIM1 on Th2 cells.²⁶ We therefore hypothesize that TIMD4 is an important immune regulator in TO. Overexpression of TIMD4 in TO may induce innate immune cells such as macrophages and dendritic cells to phagocytose apoptotic T cells to maintain immune balance in the context of increased apoptotic load during the active phase of the disease. Increased TIMD4 expression may also explain the Th2 dominance observed in the late phase of TO.²⁷ The hypothesis that TIMD4 overexpression is a precursor for switching off active orbital inflammation in TO is now worth exploring. 

We report here significantly higher levels of α-defensins (DEFA1 and DEFA3) in active TO. Thyroid orbitopathy may share these biomarkers with other autoimmune disorders. There is considerable evidence that α-defensins are part of the innate immunity derived from neutrophils and monocytes.^28–30 α-defensins are involved in the inflammatory response as a chemotactic factor for neutrophils, monocytes, and T cells and promote adaptive immune responses. Conversely, α-defensins also reduces inflammatory injury through inhibiting lipopolysaccharide-mediated responses, and inhibiting proinflammatory cytokine secretion by macrophages.²⁹ The role of DEFA1-3 in the pathogenesis of autoimmune diseases such as rheumatoid arthritis (RA), SLE, Behçet's disease, and lung diseases has been increasingly recognized.^{28,29,31–33} 

The α-defensin-1 promotes inflammation through induction of proinflammatory cytokines IL-6, IL-8, and matrix metalloproteinases, possibly through regulation of JNK and/or ERK MAP kinases and NF-κβ pathways.³¹ Furthermore, α-defensins induce synthesis of mucin and chemokines including IL-8, TNF-β, and VEGF in lung epithelial cells and lung fibroblasts cultures.²⁸ Even though α-defensins are mainly expressed by neutrophils, other studies suggest that T lymphocytes, NK cells, and monocytes also produce DEFA1-3. DEFA1-3 are known to be significantly higher in T cells from peripheral blood and blister fluid in Stevens-Johnson syndrome/toxic epidermal necrolysis patients, where cytotoxic T cells are the main effector cells.³⁰ 

Our study confirms findings of previous microarray and cell culture studies that increased adipogenesis occur in active TO.^6,34,35 A group of genes involved in fatty acid and cholesterol synthesis was upregulated in our study including SREBF1, SCD, FADS1, and SCDP1. Gene sets enriched in lipid metabolism, including biosynthesis of unsaturated fatty acids, fatty acyl-CoA synthesis, triacylglycerol and ketone body metabolism, glycolysis pathway, and integration of energy metabolism, were upregulated. In addition, pathway analysis showed overrepresentation of pathways linked to adipogenesis including PPAR activation, insulin signaling pathway, unsaturated fatty acid synthesis, and glycerolipid metabolism. 

The overexpression of SCD in TO orbital fat tissue compared with normal controls was independently validated by Lantz et al.,⁷ who further linked overexpression of immediate early genes in smokers with severe active TO.³⁶ Upregulation of several adipocyte regulatory genes including PPARγ, adiponectin, and leptin support increased de novo adipogenesis and/or fat cell expansion in TO.⁶ In addition, adipogenesis may be regulated through Wnt-signaling with dysregulation of secreted frizzled-related protein-1, Wnt5a, sFRPs, and DKK.^6,37 Our study identified frizzled family receptor 9 (FZD9), a receptor for Wnt signaling protein, which was modestly downregulated in active compared to inactive TO. FZD9 was activated by Wnt-2 and functions in Wnt/β-Catenin signaling in mice.³⁸ However, a pathogenic role for FZD9 in TO remains unclear as yet. 

Multiple cell-mediated immune genes were overexpressed in TO confirming the importance of cellular immunity in the pathogenesis of TO. CD247, CD3D, CD8A, CST7, GZMA, CAMP, SLAMF6, EOMES, LTB, and CCL5 were upregulated in active TO compared with inactive TO. CD3d and CD247 are components of the T-cell receptor complex, while Cd8a is a surface molecule on cytotoxic T cells. CST7 is a member of cysteine protease inhibitors expressed by T cells, natural killer cells, monocytes, and mast cells.³⁹ CST7 could be important in regulating serine protease activity in IL-12 activated natural killer cells by inhibiting cathepsin C.⁴⁰ Granzyme A (GZMA) is a serine protease found in cytotoxic T cells and NK cells. There is mounting evidence that granzyme A induces caspase-independent apoptosis in targeted cells, as well as inducing IL-8, IL-6, and IL-1β production and facilitates lymphocyte migration through the extracellular matrix.⁴¹ SLAMF6 encodes a transmembrane protein of CD2 subfamily on NK cells, T and B lymphocytes, and functions as a costimulatory molecule, which primes T cells to produce Th1 cytokines.⁴² CAMP relates in function to defensin in antimicrobial activity, cell chemotaxis, immune mediation, and inflammation response regulation.²⁹ PTPRCAP is a positive regulator of tyrosine phosphatase PTPRC (CD45), a key regulator of T and B lymphocytes activation.⁴³ EOMES is a transcription factor crucial for mesoderm and neural development and a key regulator for differentiation of CD8+ T cells.⁴⁴ LTB, also known as TNF superfamily, member 3 is proinflammatory and is involved in normal development of lymphoid tissue.⁴⁵ The upregulation of LTB supports the notion of a role of TNF-antagonist in the treatment of active TO.⁴⁶ GPX3 functions in detoxification of reactive oxygen species and hydrogen peroxide that involves selenium binding at its active site (http://www.ncbi.nlm.nih.gov/gene/2878; provided in the public domain by the National Center for Biotechnology Information, Bethesda, MD, USA). The upregulation of GPX3 in active TO may explain the observations that selenium supplementation helps with soft tissue signs in active TO, and the relative selenium deficiency in TO cases.^47,48 Of note, targets for novel and experimental therapies for TO including CD-20, IL-1, IL-6, IL-8, and IGF-1R were not differentially expressed in the orbital adipose tissue that we studied, and TSHR was low ranking and marginally downregulated in active TO compared to inactive TO (data not shown). 

CCL5, also known as RANTES, is an important chemo-attractant for blood monocytes, T-helper cells, and eosinophils. Immunoglobulins from Graves' disease patients induced increased expression of RANTES in TO orbital fibroblast.⁴⁹ The GSEA and pathway analysis are in keeping with translational studies highlighting the Th1 cytokine expression profile, T-cell recruitment, and activation by orbital fibroblasts, dysregulation of MAPK signaling, infiltration of lymphocytes, macrophages, and monocytes in TO affected orbital tissues on histology.^5,49–53 

Our study revealed epigenetic factors may play a role in pathogenesis of TO. Active TO displayed gene set enriched with epigenetic signatures of acute myeloid leukemia (AML). Hypermethylation of DNA leads to silencing of genes, conversely aberrant hypomethylation of cluster 5 epigenetic markers in AML relates to dysregulation of immunity-related pathways involving immune deficiency signaling, cytotoxic T-cell mediated apoptosis, and T-cell receptor signaling.⁵⁴ Active TO also showed downregulation of genes with unmethylated histone 3 in high CpG density promoters (HCP). Methylation at K4 (H3K4me3) and K27 (H3K27me33) on histone 3, known as bivalent status on HCP is associated with cell lineage commitment, genes related to adipogenesis (PPARγ) often remains bivalent in differentiating embryonic fibroblast.⁵⁵ 

Gene sets enriched with T cell-mediated immune response include T-cell receptor signaling via COOH-terminal Srk kinase (CSK), Lck and Fyn tyrosine kinases, TOB1 nuclear transducer together modulate T-cell activation, and proliferation. Gene sets involved in IL-12 and NO2 dependent IL-12 signaling pathways were upregulated, which induce Th1 cell differentiation and activation of both T and NK cells via JAK/STAT transcription factor signaling. Activation of T-helper cell, cytotoxic T cell, caspase cascade, and telomere shortening mediated immune-related apoptosis in target cells. Further gene sets enriched in cell-to-cell adhesion signaling, stathmin, and E-cadherin signaling mediating intercellular communication may result in changes in cytoskeleton and facilitate cell mobility, migration, and proliferation. 

One potential limitation of the interpretation of the microarray data is that prior immunosuppressive treatment was not factored into the analysis. Despite more active TO cases had intravenous methylprednisolone within 2 months prior to surgery compared to inactive TO and normal controls, active TO cases still had signs of orbital inflammation at the time of surgery. The microarray results still detected significantly upregulation of genes involved in immune and inflammatory responses and showed positive correlation of inflammatory index score in vast majority of the top ranked genes expressions tested. In addition, intravenous dexamethasone was routinely used at induction of general anesthesia for orbital surgery in both active TO and inactive TO and normal controls, and the prevalence of intravenous dexamethasone use between the groups of patients was statistically insignificant; hence, the differential RNA expression observed is unlikely due to confounding steroid effects. 

In conclusion, TIMD4, DEFA1, DEFA1B, and DEFA 3 were overexpressed in active TO compared with inactive TO, suggesting a pathogenic role of the innate immune response in TO. Active TO was marked by upregulation of genes involved in cell-mediated, innate, and inflammatory responses with concurrent enhancement of orbital adipogenesis. Epigenetic mechanisms may play a role in TO. The study gives new insights into the underlying complex molecular mechanisms in TO, and provides novel insights into candidate molecules and pathways, which can be explored to develop alternative treatment strategies for the treatment of TO, which carries substantial visual morbidity. 

Acknowledgments

Supported by Ophthalmic Research Institute of Australia (ORIA) new investigator grant. The authors alone are responsible for the content and writing of the paper. 

Disclosure: J.J. Khong, None; L.Y. Wang, None; G.K. Smyth, None; A.A. McNab, None; T.G. Hardy, None; D. Selva, None; B. Llamas, None; C.-H. Jung, None; S. Sharma, None; K.P. Burdon, None; P.R. Ebeling, None; J.E. Craig, None