Gene expression analysis of the cerebellar cortex in essential tremor

Journal Pre-proof Gene Expression Analysis of the Cerebellar Cortex in Essential Tremor Regina T. Martuscello, Chloe¨ A. Kerridge, Debotri Chatterjee, Whitney G. Hartstone, Sheng-Han Kuo, Peter A. Sims, Elan D. Louis, Phyllis L. Faust

PII:

S0304-3940(19)30643-3

DOI:

https://doi.org/10.1016/j.neulet.2019.134540

Reference:

NSL 134540

To appear in:

Neuroscience Letters

Received Date:

2 July 2019

Revised Date:

28 September 2019

Accepted Date:

8 October 2019

Please cite this article as: Martuscello RT, Kerridge CA, Chatterjee D, Hartstone WG, Kuo S-Han, Sims PA, Louis ED, Faust PL, Gene Expression Analysis of the Cerebellar Cortex in Essential Tremor, Neuroscience Letters (2019), doi: https://doi.org/10.1016/j.neulet.2019.134540

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier.

Gene Expression Analysis of the Cerebellar Cortex in Essential Tremor Regina T. Martuscelloa, Chloë A. Kerridgea, Debotri Chatterjeea, Whitney G. Hartstonea, ShengHan Kuoc, Peter A. Simsd,e,f, Elan D. Louisg,h,i, Phyllis L. Fausta a

ur na

Corresponding author: Phyllis L. Faust, MD, PhD [email protected] Tel: (212) 305-7345

lP

re

-p

ro

of

Department of Pathology and Cell Biology, Columbia University Medical Center and the New York Presbyterian Hospital, 630 W 168th Street, P&S 15-405, New York, NY, USA. b College of Physicians and Surgeons, Columbia University Medical Center and the New York Presbyterian Hospital, 630 W 168th Street, New York, NY, USA. c Department of Neurology, College of Physicians and Surgeons, Columbia University, 630 W 168th Street, BB302 New York, NY, USA. d Department of Systems Biology, Columbia University Medical Center, 3960 Broadway, RM208, New York, NY, USA. e Sulzberger Columbia Genome Center, Columbia University Medical Center, 1150 St. Nicholas Ave, New York, NY, USA. f Department of Biochemistry & Molecular Biophysics, Columbia University Medical Center, 701 W 168th Street, New York, NY, USA. g Department of Neurology, Yale School of Medicine, Yale University, 15 York Street, New Haven, CT, USA. h Department of Chronic Disease Epidemiology, Yale School of Public Health, 15 York Street, Yale University, New Haven, CT, USA. i Center for Neuroepidemiology and Clinical Neurological Research, Yale School of Medicine, Yale University, 15 York Street, New Haven, CT, USA.

Jo

Regina T. Martuscello, PhD: [email protected] Chloë A. Kerridge, MSc: [email protected] Debotri Chatterjee, BA: [email protected] Whitney G. Hartstone, BA: [email protected] Sheng-Han Kuo, MD: [email protected] Peter A. Sims, PhD: [email protected] Elan D. Louis, MD, MSc: [email protected]

Highlights

    

Essential tremor (ET) is one of the most common neurological diseases The underlying molecular mechanisms are not clear Our goal was to explore the underlying molecular source of ET disease pathogenesis We performed a transcriptomic analysis by direct sequencing of RNA from cerebellum The analyses draw attention to several genes and regulatory pathways

Abstract

of

Essential tremor (ET) is one of the most common neurological diseases, with a central feature of an 8-12Hz kinetic tremor. While previous postmortem studies have identified a cluster of

ro

morphological changes in the ET cerebellum centered in/around the Purkinje cell (PC)

population, including a loss of PCs in some studies, the underlying molecular mechanisms for these changes are not clear. As genomic studies of ET patients have yet to identify major genetic

-p

contributors and animal models that fully recapitulate the human disease do not yet exist, the study of human tissue is currently the most applicable method to gain a mechanistic insight into

re

ET disease pathogenesis. To begin exploration of an underlying molecular source of ET disease pathogenesis, we have performed the first transcriptomic analysis by direct sequencing of RNA

lP

from frozen cerebellar cortex tissue in 33 ET patients compared to 21 normal controls. Principal component analysis showed a heterogenous distribution of the expression data in ET patients that only partially overlapped with control patients. Differential expression analysis identified 231

ur na

differentially expressed gene transcripts (‘top gene hits’), a subset of which has defined expression profiles in the cerebellum across neuronal and glial cell types but a largely unknown relationship to cerebellar function and/or ET pathogenesis. Gene set enrichment analysis (GSEA) identified dysregulated pathways of interest and stratified dysregulation among ET cases. By GSEA and mining curated databases, we compiled major categories of dysregulated processes

Jo

and clustered string networks of known interacting proteins. Here we demonstrate that these ‘top gene hits’ contribute to regulation of four main biological processes, which are 1) axon guidance, 2) microtubule motor activity, 3) endoplasmic reticulum (ER) to Golgi transport and 4) calcium signaling/synaptic transmission. The results of our transcriptomic analysis suggest there is a range of different processes involved among ET cases, and draws attention to a particular set of genes and regulatory pathways that provide an initial platform to further explore the underlying biology of ET.

Key Words: Essential Tremor, Cerebellum, RNA-seq, Purkinje Cell, GSEA, String networks

Introduction

of

Essential tremor (ET) is a common neurological disease whose primary clinical feature is an 812Hz kinetic tremor of the arms; the disease affects approximately 4% of adults over 40 and an

ro

estimated excess of 20% of adults over 90 [1]. Patients diagnosed with ET are predominantly adults, but tremor can occur earlier in life as well, with younger patients predominantly having

-p

the familial form of the disease [2]. ET is progressive in nature, with patients experiencing worsening tremor [3, 4]. Diagnosing ET can be challenging [5, 6], as tremor is a common feature

re

of other neurological diseases, there are currently no disease biomarkers, and tremor, in itself, is commonly noted in the normal elderly [7]. Traditionally, little has been known about the

lP

biological basis for ET and few autopsies were performed. Scientific investigation places the cerebellum with its oscillatory network comprised of the cerebello-thalamo-cortical loops at the center of ET research [8, 9], but a definitive etiology for tremor still does not exist. To advance

ur na

understanding of the underlying patho-mechanisms of ET, in 2003, we established the Essential Tremor Centralized Brain Repository (ETCBR) [10]. To date, this has resulted in the prospective harvest and postmortem study of approximately 200 ET brains [11]. Through this mechanism, over the 15 years, we have carefully catalogued a cluster of morphological changes in the ET cerebellum centered in/around the Purkinje cell (PC) and

Jo

neighboring neuronal populations [12, 13]. We have detailed specific alterations in PC axonal and dendritic morphology in ET patients compared to age-matched controls, including increased numbers of torpedoes (ovoid axonal swellings) [14], thickened axonal profiles, branching axons and recurrent axonal collaterals [15]; increased numbers of dendritic swellings [16], regressive changes in the PC dendritic architecture (reduced dendritic complexity and dendritic spine density) [17], decreased climbing fiber synaptic density on PC dendritic shafts and increased extension of climbing fibers into the parallel fiber territory in the molecular layer of the dendritic

arbor [18]. A significant loss of PCs has been identified in several although not all studies [12]. A number of other changes, including increases in the number of heterotopic PCs [19, 20], “empty baskets” (a marker of PC loss) [21] and hypertrophic basket cell axonal processes [22] have been noted in ET patients versus controls. Cerebellar expression of the excitatory amino acid transporter 2 (EAAT2), a major glial glutamate transporter, was differentially altered with decreased expression in cerebellar cortex and increased expression in the dentate nucleus of ET patients versus controls [23]. These changes represent an array of pathological deviations from

of

the norm, and along with heterogeneity in clinical findings, suggest that ET may very well be a family of related diseases rather than a single clinical-pathological entity.

Something must set this degenerative cascade of events into motion in the ET cerebellum.

ro

Indeed, one can hypothesize the presence of an as-yet unidentified molecular event or events that results in a breakdown in PC house-keeping functions and cellular physiology, with the observed

-p

cellular changes representing a series of responsive downstream markers of neuronal duress. To shed light on the underlying molecular basis, we performed a high-throughput RNA-sequencing

re

of whole cerebellar cortex lysates from 33 ET patients compared to 21 normal age-matched controls. Despite a heterogeneous profile of cerebellar gene expression changes among ET

lP

patients, differential expression analysis, curated database mining and gene set enrichment analysis (GSEA) have identified four main biological processes that provide new leads into the

Methods

ur na

underlying molecular basis of cellular dysfunction in ET.

Study design and human sample collection ET brains obtained from the Essential Tremor Centralized Brain Repository (ETCBR), a longstanding joint collaboration between investigators at Yale and Columbia Universities [10].

Jo

This centralized repository provides brains of patients with ET in the United States. ET diagnoses were carefully confirmed by a senior movement disorder neurologist (E.D.L.) using 3 sequential methods [11, 15]. Briefly, the clinical diagnosis of ET was initially assigned by treating neurologists, and secondly, confirmed by E.D.L. using questionnaires, review of medical records and review of Archimedes spirals. Third, a detailed, videotaped, neurological examination was performed. Combined with the questionnaire data, the final diagnosis of each ET case was re-examined, using previously published diagnostic criteria, which have been

shown to be both reliable and valid [24]. None of the ET cases had a history of traumatic brain injury, exposure to medications known to cause cerebellar damage, or heavy ethanol use. Every 6 months a follow-up semi-structured telephone evaluation was performed and hand-drawn spirals were collected until the death of the patient [21]. Eight control brains were from the New York Brain Bank. Individuals were prospectively followed at the Alzheimer’s Disease Research Center or the Washington Heights Inwood Columbia Aging Project at Columbia University [1]. During serial neurological examinations, these individuals were clinically free of ET and other neurodegenerative disorders, including

of

Alzheimer’s disease, Parkinson’s disease, or progressive supranuclear palsy. Seven control

brains obtained from the National Institutes of Health NeuroBioBank (three from University of

ro

Maryland, Baltimore, MD and four from University of Miami, Miami, FL); five control brains were obtained from the University of Kentucky Alzheimer’s Disease Center Brain Bank

re

respective university or institutional ethics boards.

-p

(Lexington, KY). All study subjects agreed and signed informed consent forms approved by the

Tissue Processing

lP

All ET and 8/21 control samples had a complete neuropathological assessment at the New York Brain Bank as previously described [21]. Briefly, all brains had standardized measurements of brain weight (grams) and post-mortem interval to freezing. Paraffin blocks from all brain regions

ur na

stained with standard pathology series (Luxol fast blue/hematoxylin and eosin [LH&E], Bielschowsky, -synuclein, -amyloid, phosphorylated tau [AT8]) and examined microscopically, including Braak and Braak Alzheimer’s disease staging for neurofibrillary tangles, and Consortium to Establish a Registry for Alzheimer’s disease (CERAD) ratings for neuritic plaques. Cerebellar control tissue (13/21) from outside brain banks were harvested from

Jo

the same region as at the New York Brain Bank. For most cases, brain weight (grams) was available, and in some instances, data was also available on postmortem interval and Braak Alzheimer’s disease staging. All frozen cerebellar tissue blocks were obtained from a parasagittal slice located 1-1.5 cm from the cerebellar midline, containing anterior and posterior quadrangulate lobules, a region known to be involved in motor control. A standard cerebellar tissue block from the same region was paraffin embedded and PCs were quantified in a 7-m-thick, LH&E-stained section by counting and averaging the number of PCs

across fifteen 100x non-overlapping microscopic fields. Similarly, torpedoes (i.e. fusiform swellings of the proximal portion of the PC axon) were counted in one entire LH&E section, and the resulting total torpedo counts were normalized to PC layer length (i.e. divided by PC layer length) to account for any potential variations in amount of cerebellar cortex in the tissue block.

of

Sample Selection Samples selected for the study underwent additional clinical and pathologic screening prior to inclusion. Demographic information for all selected samples is provided in Supplemental Table

ro

1. ET and control brains were chosen based on a lack of widespread marked hypoxic-ischemic damage, concurrent Alzheimer's-type changes that would meet criteria for high likelihood of

-p

Alzheimer's disease, or evidence of other neurodegenerative disease pathology (e.g., progressive supranuclear palsy, corticobasal degeneration, Lewy bodies, traumatic encephalopathy). ET

re

patients that underwent deep brain stimulation surgery were excluded. A senior neuropathologist (P.L.F) examined LH&E stained cerebellar slides to determine good tissue morphology,

lP

including absence of autolytic changes in the granule cell layer. Samples were then age-matched. There were 55 control samples overall and after above assessments ended with 22 sample library preparations, 21 samples sequenced (1 library not usable for sequencing) and 20 samples utilized

ur na

in the overall analysis (1 under-sequenced sample removed as outlier). There were 97 ET samples overall and after above assessments ended with 33 sample library preps, 33 samples sequenced and 32 samples utilized in the overall analysis (1 under-sequenced sample removed as outlier).

Jo

RNA Extraction & Quality Control

All samples had total RNA extraction and RNA integrity numbers (RINs) assessed prior to sequencing. Briefly, a 1cm2 block containing only cerebellar cortex was sampled and three100m thick cryostat sections were placed in a 2 mL round bottom RNase free tube that was precooled on dry ice. RNA extraction was performed at the Columbia University Molecular Pathology Core utilizing a QIASymphony instrument. RINs were evaluated using an Agilent Bioanalyzer 2100 Nano Chip. Among sixty total samples (36 ET & 24 control), RINs excluded

ranged from 4.2-6.3 (5 samples [3 ET, 2 control]; average 5.5) and RINS accepted ranged from 7.8-10 (55 samples; average 9.3, Supplemental Table 1).

Library preparation & RNA sequencing All samples were prepared simultaneously to avoid batched library preparation bias using modifications from previously published protocols [25]. Briefly, mRNA was isolated from 450ng of total RNA by incubation with streptavidin-coated C1 Dynabeads (ThermoFisher

of

CAT#65001). Dynabeads were washed (20mM Tris pH 8.0, 50mM NaCl, 0.1% Tween20/0.5mL per 20uL of beads), placed in hybridization buffer (20mM Tris pH 8.0, 1M NaCl,

0.1% Tween-20/20uL per 20uL of beads) and coated with 25uM biotinylated LNA-Oligo(dT)

ro

(Exiqon CAT#300100-03/4uL per 20uL of beads). LNA-Oligo(dT)-coated Dynabeads were incubated with each sample at 1:1 volume with RNase inhibitor added (SUPERaseIN

-p

[ThermoFisher CAT#AM2694], 0.50 l per sample). Libraries were constructed using SMARTer Stranded RNA-Seq kit (Clontech Laboratories CAT#634862) according to the

re

manufacturer’s instructions. mRNA was fragmented for 4 minutes at 94°C and cDNA libraries were purified using AMPure XP beads (Beckman Coulter CAT#A63880). Libraries were

lP

amplified via 9-cycles of PCR and each sample indexed with unique barcodes for highthroughput sequencing. RNA-seq library was purified again using AMPure XP beads. The cDNA quality in the sequencing libraries was determined using an Agilent Bioanalyzer 2100

ur na

High Sensitivity DNA assay and concentration determined in Qubit assay (ThermoFisher CAT#Q32851). The libraries were pooled (at equal concentrations) and 1x75 base-pair singleend sequenced in quadruplicate on an Illumina NextSeq500 using 75-cycle High Output Kits (~30M reads per sample) at the Columbia University Genome Center. Raw sequencing data and

Jo

DEseq2 analysis can be found at the NIH Gene Expression Omnibus (GEO): GSE134878.

RNA sequencing Analysis We aligned the RNA-seq reads to the human genome and transcriptome annotation (hg19, UCSC annotation from Illumina iGenomes) using STAR aligner (STAR v2.4.2a) and quantified the reads that aligned uniquely to the transcriptome using the feature Counts program in the subread package (subread v1.4.6). Mean alignment rate was 83% (Supplemental Table 1). We identified differentially expressed genes using DESeq2 [26]. We conducted principal component

analysis (PCA) of the RNA-seq data set by first normalizing the counts for each sample using size factors obtained from DESeq2, log-transforming the normalized counts, and computing the principal components of the resulting matrix using the scikit-learn module in Python [27].

Curated Databases We selected 231 ‘top’ differentially expressed genes (DEseq2 p < 0.025, false discovery rate < 0.25) and examined protein expression in The Human Protein Atlas [28], and gene pathway

of

analyses using the STRING network database [29], WebGestalt [30], Reactome [31] and The PANTHER (Protein ANalysis THrough Evolutionary Relationships) Classification System [32].

ro

Gene Set Enrichment Analysis (GSEA)

We used Gene Set Enrichment Analysis (GSEA – Java implementation from the Broad Institute,

-p

http://software.broadinstitute.org/gsea/index.jsp) to assess the enrichment of gene ontologies in the principal components [33, 34]. We first computed the Spearman’s correlation coefficient

re

between the relative abundance of each gene across all samples and the first two principal components (PC1, PC2). For each principal component, we ranked all genes by their Spearman’s

lP

correlation coefficient and used this as input for GSEA (pre-ranked /classic mode). To obtain the normalized enrichment score for gene ontologies, we used the C5 curated gene sets from the Molecular Signatures Database (GO Biological Process, GO Cellular Component, and GO

ur na

Molecular Function). For single-sample GSEA (ssGSEA), we first log-normalized by the total number of reads in each sample and then computed z-scores for each gene. We used the z-scores to rank genes in each sample, and these ranked lists were used as input to GSEA (preranked/classic mode) with the C5 curated gene sets from the Molecular Signatures Database.

Jo

Results

Heterogeneity of ET Cerebellar Transcripts Gene transcripts were analyzed across 33 ET samples and 21 controls (NIH-GEO: GSE134878). First, we performed a principal component analysis (PCA) of the gene expression profiles across all 54 samples to determine whether there was diagnosis-specific clustering or outliers in the sample set (Supplemental Figure 1A; blue dots=control; red dots=ET). Two outlier samples were identified and removed from subsequent analyses, as they were found to be

under-sequenced, having the lowest number of reads per gene of any sample (Supplemental Figure 1A; black circle, 1 control [50% of genes], 1 ET [40% of genes]). PCA was then reassessed without outliers (Figure 1A). Tight clustering according to clinical diagnosis was not observed. However, there was only partial overlap between ET and control samples, with ET samples having a significantly higher PC2 value compared to controls (Figure 1B; ET mean = 9.13, median = 7.18; control mean = -14.61, median = -22.77; p=0.020 two-tailed T-Test). Overall, the relative expression patterns across the ET samples show variability among the ET

of

patient population tested.

To investigate potential diagnostic differences between ET and control samples, we performed

ro

differential gene expression analysis comparing 32 ET samples versus 20 controls with DEseq2 (NIH-GEO: GSE134878). The relative log2fold changes observed between control and ET

-p

samples were modest (Supplemental Table 2, Supplemental Table 3). We used a FDR cutoff value of 25% (<0.25) to expand the gene pool for study in the context of an exploratory dataset

re

to uncover potential hypotheses for future research validation [34]. This resulted in 231 top gene

lP

hits.

To understand the scientific context in which these genes relate to the study of ET and the cerebellum, we thoroughly vetted the top gene hits using The Human Protein Atlas [35] and

ur na

PubMed. Using The Human Protein Atlas, we compiled immunohistochemical staining localization within the cerebellar cortex for each of the 231 top hits [columns 7-9], as well as previously published RPKMs (reads per kilobase million) from the cerebellum and the cerebellar cortex [columns 10 & 11] (Supplemental Table 3). Using PubMed, we compiled a publication list for genes relevant in the study of ET by combining search terms located in the title, abstract,

Jo

or keywords of a publication for each gene with Purkinje cell, GABA (gamma-Aminobutyric acid; PCs are GABAergic), cerebellum and brain [columns 12-16] (Supplemental Table 3). Overall, there is a paucity of scientific evidence available as to how these genes may or may not regulate ET, the cerebellum and the human brain, as indicated by the large red block of 0-1 available publications in these areas [columns 12-16] (Supplemental Table 3). There is, however, significant information on the basic regulatory biological role a subset of these genes

play in the human body, the brain and in other neurodegenerative diseases (e.g., UBQLN1, KLC2, DNN3, S100B, SYNJ1, DYNC1H1, TTBK2, UBTF, CAPN2).

Dysregulated Clusters of String Networks To further investigate the potential role of the top gene hits, we used online curated databases to categorize the regulatory interactive functions of these genes in the human body, as well as in disease states. We employed a combinatorial approach utilizing Web-based Gene Set Analysis

of

Toolkit (Web-Gestalt), Reactome and Panther databases to reveal biological insights from the functional enrichment analysis software and identify significantly enriched pathways (Table 1). This classification process delineates multiple categories of biological functions relevant to the

ro

nervous system regulated by these genes, with localization across the cell body, axon and dendrite, and significant regulation of several homeostatic cellular functions. Overall, we

-p

identified four main categories of dysregulation that numerous processes fall into: (1) axonal guidance, (2) microtubule motor activity, (3) ER-Golgi transport, and (4) calcium

re

signaling/synaptic transmission.

lP

To visualize how these biological processes and regulatory signaling pathways relate to one another, we applied the String Network database for a system-wide grasp of functional interactions between expressed proteins [29]. Figure 2A shows 230/231 (LINC00641 transcript

ur na

did not map; non-coding RNA) expressed proteins clustered by the top 7 k-mean groups [axon guidance (red), kinesin/microtubule (dark blue), ribosomal (yellow), p53/cell cycle (light blue), ER-Golgi transport (green), calcium signaling (purple); all others in tan]. Clustering by k-means is a method of vector quantization, whereby the clustering aims to screen n number of observations into k groupings in which each observation belongs with its nearest neighbor.

Jo

Notably, there are numerous proteins that have no known or hypothesized interactions with any of the other well-studied networks, confirming the current lack of scientific knowledge for many of our top gene hits. We then focused on four identified pathways that show potential relevance in the field of ET research, which are separately clustered in Figure 2B to show their interactive significance [axon guidance (red), kinesins (blue), ER-Golgi transport (green), calcium signaling (purple)]. This interactive clustering and relationship among the relevant pathways indicate that dysregulation in just one of these proteins can have significant effects on numerous other

pathways. Taken together, it is increasingly likely that multiple dysregulated processes play a role in the etiology of ET, further supporting the hypothesis that ET is a family of related diseases rather than a single entity.

Gene Set Enrichment Analysis (GSEA) The splay shown across the ET and control samples in the PCA (Figure 1A) suggest a level of heterogeneity in the sample population tested. To investigate biological sources of variability

of

across these samples and better understand the potential biological relevance of these data, we employed a two-pronged Gene Set Enrichment Analysis (GSEA) to first analyze gene expression in all samples based on the principal components and then again as independent variables. This

ro

analytic approach enabled us to investigate (1) clusters of dysregulation by diagnosis and (2) clusters of dysregulation among the ET patients. Table 2 displays positively and negatively

-p

correlated pathways that are enriched in principal component 1 (PC1) and principal component 2 (PC2) for the entire dataset, and is broken down into Biological Process, Cellular Component

re

and Molecular Function for each. Notably, all four of our main categories of dysregulation are present in the GSEA, divided across PC1 and PC2. ET samples are enriched in positive PC2

lP

space compared to controls, and all samples distribute across PC1 (Figure 1B). Gene ontology (GO) pathways for positive PC2 correlations are enriched for axonal, dendritic, synaptic and post-translational alterations. GO pathways for positive PC1 correlations overlap with some

ur na

positive PC2-correlated pathways, but are enriched for ER-Golgi transport vesicles and protein localization. GO pathways for various molecular level-binding activities (selective, non-covalent, stoichiometric interactions) are enriched in both PC1 and PC2. Given that a significant portion of ET cases fall in positive PC2 space and enriched pathways for PC2 involve axonal and microtubule/kinase/motor activity we examined potential correlations between pathology metrics

Jo

for axonal changes with PC2 values. We found a significant correlation between PC2 and torpedo counts per Purkinje cell (Supplemental Figure 2A; Spearman p=0.0049 r=0.3845). As a whole, this dataset analysis shows a heterogeneous mix of dysregulated processes that divides samples into four main biological categories: (1) axon guidance, (2) microtubule/kinase motor activity, (3) ER-Golgi transport and (4) calcium signaling/synaptic transmission.

The GSEA enrichment scores represent the phenotypic differences across the entire dataset. To further examine the heterogeneity seen in our ET samples, we employed a single-sample GSEA (ssGSEA), calculating separate normalized gene set enrichment scores (NES) for each sample’s gene expression profile (Figure 3A). The NES characterizes the activity level of each GO process, whereby the members in each gene set are coordinately up- (red) or down- (green) regulated. Hierarchical clustering of the NES for each GO pathway (Table 2) demonstrated statistically significant variation and two main groupings of ET patients in the heatmap (Kruskal-

of

Wallis H test p=0.0001), whereby some patients have significantly more gene expression within the presented pathways of dysregulation than others (Figure 3A). The spearman correlation

matrix, which is calculated coefficients (r values) between each ET sample, shows significant

ro

positive (blue) and negative (red) correlations between these two ET patient groupings (Figure 3B), with some cases straddling each group (white = no correlation). Taken together, this

-p

demonstrates degrees of biological dysfunction occurring within the ET patient population examined and this spectrum of dysfunction supports our claim that ET is a family of related

re

disorders.

lP

Discussion

ET is one of the most common neurological diseases in adults, with a high degree of clinical heterogeneity that has made gold standard diagnostics challenging. To date, genetic studies have

ur na

focused on ET families, as many affected individuals have familial disease history [36]. Even with this strong hereditary component, only a handful of pathogenic risk variants have been reported, and none has shown significant pathogenicity in vivo or reproducibility in patient populations tested. It should be noted that familial ET was found to be pathologically equivalent to sporadic cases of ET when examined in postmortem cerebellum [37]. Genome-wide

Jo

association studies identified LINGO1 and SLC1A2 as potential ET susceptibility genes in Icelandic & North American populations and European populations, respectively [38]. Notably, we (and another [39]) did not see any significant differential expression in these genes in our dataset. Exome sequencing of ET families identified potential pathogenic variants of FUS, HTRA2, TENM4 and NOS3 as possibly causative genes [40]. Genome-wide scans and linkage analyses of large affected families have identified 3 genetic loci (ETM 1-3) of chromosomal abnormality [39, 41]. Whole genome sequencing (WGS) of early onset families also identified

NOS3, as well as KCNS2 (Kv9.2), HAPLN4 (BRAL2) and USP46 as gene variants predicted to affect function [42]. KCNS2 (Kv9.2) has shown some disease recapitulation in a Drosophilia hKv9.2 transgenic model [43]. Most recently, another WGS of ET families identified additional candidate genes, including a deleterious damaging variant in CACNA1G (Cav3.1), which encodes the pore forming subunit of the T-type calcium channels; SLIT3, which encodes an axon guidance molecule; and in three families phenolyzer prioritized genes that were associated with hereditary neuropathies (Family A – KARS, Family B – KIF5A, Family F – NTRK1) [44].

of

Interestingly, axonal guidance disruption and dysregulation of calcium signaling/synaptic transmission were two of the four identified pathways from our investigation into the top gene hits, and KIF5A was found to be a statistically significant finding in our experimental groups,

ro

along with numerous other kinesin genes. While our patient population has not been tested for the same genetic variants as these ET families, the physiological roles of our genes of interest are

-p

in line with some familial studies, as well as the spectrum of pathological changes seen in the ET

Axon Guidance, Microtubules and Kinesins

re

patient cerebellum.

lP

Our pathologic studies in ET and emerging evidence in other neurodegenerative diseases support the dysregulation of pathways involved in axonal maintenance. Among the broad range of structural and degenerative changes we have described in the ET cerebellum, a plethora of those

ur na

changes are involved in PC axonal morphology, including an increased number of thickened axonal profiles, torpedoes, axonal branching, axonal recurrent collaterals and terminal axonal sprouting [45]. Electron microscopy studies have shown that axonal swellings in ET patients contain an excess of disorganized neurofilaments and organelles (mitochondria and ER) [46], which accumulate as a result of damaged axonal transport [15]. Trafficking along the axon

Jo

requires stable microtubule dynamics via dynein-dynactin (retrograde) and kinesin (anterograde) activity; the need to maintain tight regulation of neuronal microtubule dynamics are now emerging, as defects have been shown in a several neurodegenerative disorders [47]. We saw numerous genes in our top hits that regulate these important processes in the nervous system. Most notably, all 13 genes involved in kinesin anterograde transport were upregulated in our ET cases compared to controls (KIF5A, KIF5B, KIF2A, KLC1, KLC2, DYNC1H1, RAB3GAP1, RAB6B, TTBK2, SPTBN5, SPTAN1, COPB2, ARCN1).

ER-Golgi Transport Additional control of microtubules and microtubule-associated proteins comes from the regulatory function of the Golgi apparatus-complex that has a central role in transport, processing and sorting of proteins. Fragmentation of the neuronal Golgi apparatus has been reported in multiple neurodegenerative diseases (Amyotrophic lateral sclerosis, Spinocerebellar ataxia 2, Alzheimer’s Disease) as an early and irreversible lesion resulting from mechanisms that

of

revolve around microtubule stabilization and stress [48]. It has also been shown that Golgi fragmentation and dispersal precedes neuronal cell death, which can be triggered by excitotoxins, oxidative/nitrosative insults or ER stress [49]. The Golgi is becoming increasingly viewed as a

ro

central regulator of signaling cascades, as it spatially regulates kinases, phosphatases and

GTPases [50]. The combination of altered genes centered around the Golgi, kinesins, axonal

-p

guidance and microtubule transport in our data suggest that shifts in gene expression, and in turn

re

protein regulation, have significant effects on neurons in the human brain.

Calcium Signaling and Synaptic Transmission

lP

Most research relating to calcium channels in ET has centered around the T-type calcium channel CACNA1G (Cav3.1) in the inferior olive, and some harmaline mouse models have suggested a role in tremor [36]. It is well established that tight regulation of calcium signaling is

ur na

essential for neuronal function, as it is a major mediator of cellular excitability and a second messenger for biochemical homeostasis. PCs play a central role in cerebellar function, as they are the only efferent projection from the cerebellar cortex. Thus, as the primary output neuron in the cerebellum, excessive or insufficient calcium levels can lead to abnormal PC function, abnormal cerebellar function, abnormal firing patterns and even cell death [51].

Jo

There are numerous synapses onto PCs, including excitatory innervations from one climbing fiber per PC (CF; neuronal projections from the inferior olive) and an estimated 1014 parallel fibers (PF; axonal extensions from granule cells)-PC connections [52]. It has been shown that patients with ET have a higher percentage of CFs in the PF territory, while also showing lower CF synaptic density compared to controls [18]. CFs can generate an intermittent burst of decrementing spikes, termed complex spikes, which are slow 1-3Hz spikes that result from calcium influx into PCs, and are proposed to impose synchronous firing in PCs. It has been

hypothesized, and demonstrated in a harmaline mouse model, that increased complex spike synchrony from the inferior olive promotes tremor [53]. Other neurodegenerative diseases of the cerebellum, such as spinocerebellar ataxias, have identified several gene mutations responsible for abnormal intracellular calcium levels as well as altered cerebellar synaptic transmission and dysregulation of PC pacemaking [54]. This commonality in a genetically diverse disease has led to the postulate that altered calcium homeostasis in PCs is a possible pathological trigger in spinocerebellar ataxia. Additional support for this hypothesis comes from data showing that

of

signaling events are initiated from or lead to aberrant calcium release from the type 1 inositol 1,4,5-trisphosphate receptor (IP3R1), which is highly expressed in PCs [55]. As a whole, the

understated when discussing PCs and the cerebellum.

-p

Multiple Pathophysiologic Processes Contributes to ET

ro

biological importance of maintaining calcium signaling and synaptic transmission cannot be

It has been postulated that ET is a family of diseases, better referred to as “the essential tremors”,

re

rather than a single disease entity [56]. Looking at the spectrum of gene expression alterations identified in our patient population, and the subsequent dysregulation of essential biological

lP

processes that directly relate to ET pathology, it is probable that more than one pathophysiologic process contributes to ET, either within a given patient or across groups of patients. Cellular reactions to insult, stress or post-translational modifications can have downstream effects on

ur na

numerous areas of cellular homeostasis. This is most evident in our string networks, as we can see significant overlap in proteins from genes that regulate the four enriched pathways discussed. It is possible, that in a heterogenous disease such as ET, centered around a functionally essential cell in the cerebellum, that multiple pathways, or a spectrum of pathway dysfunction could all

Jo

lead to a similar outcome – Purkinje cell dysfunction.

Future Directions

In this study, the largest limitation is a combination of variability in the patient samples tested and observed low-fold gene expression changes that make PCR validation and other subsequent validation experiments challenging. This study utilized RNA isolated from a cerebellar cortex homogenate, which includes many different cell types along with white matter in cerebellar folia. To further elucidate potential mechanisms behind Purkinje cell changes, we are employing laser

capture microdissection of Purkinje cells in post-mortem brain for RNA-sequencing, to specifically target this vulnerable neuron in ET. Our curated database mining analysis further implicates that investigation of alterations in post-translational protein modifications are also warranted.

Conclusions Here we present a high-throughput RNA-sequencing evaluation of ET patients compared to age-

of

matched controls. We analyzed cerebellar gene transcripts using multiple methods to reveal any and all clusters of gene expression dysregulation that could explain the myriad of morphological deviations we see in ET patients. By using principal components, we first assessed the entire

ro

dataset variance to see what type of overall gene expression patterns and clusters were present. We found no significant clustering by diagnosis, although ET cases were more positively

-p

enriched in PC2, which correlated with increased torpedo and PC loss pathology. In a differential gene expression analysis comparing samples in clinical diagnostic groupings, we

re

identified specific top hit genes that distinguished ET versus control samples. With curated database mining tools, these genes were categorized into significantly enriched pathways of

lP

interest. By next employing a GSEA, we confirmed our pathways of interest and further defined a spectrum of biological, molecular and cellular processes that may lead to the downstream morphological effects seen in Purkinje cells, axons, dendrites, climbing fibers and synaptic sites.

ur na

Finally, using ssGSEA, we demonstrated the heterogeneous profile of ET samples, with significant clustering of patients into groupings representing spectrums of pathway dysregulation. These data support the hypothesis that ET is not a single entity, but rather, is a family of related disorders caused by cellular reactions from various types of biological

Jo

dysfunction, culminating in perturbations and failures in and around the Purkinje cell.

Acknowledgements The authors would like to acknowledge their funding for this project from the NIH (R01 NS088257-01A1). We would like to thank all patients that donated their brains for banking and the staff at the ETCBR at the New York Brain Bank, the NIH NeuroBioBanks and the University of Kentucky-Alzheimer’s Disease Center Brain Bank. We thank the Columbia University Genome center and Erin Bush for their assistance and expertise for this project.

References

7. 8. 9. 10. 11. 12. 13. 14. 15.

of

Jo

16.

ro

6.

-p

5.

re

4.

lP

3.

Louis, E.D., Essential tremor. The Lancet Neurology, 2005. 4(2): p. 100-110. Louis, E.D., L.N. Clark, and R. Ottman, Familial versus Sporadic Essential Tremor: What Patterns Can One Decipher in Age of Onset? Neuroepidemiology, 2015. 44(3): p. 166-72. Louis, E.D., et al., Estimating annual rate of decline: prospective, longitudinal data on arm tremor severity in two groups of essential tremor cases. J Neurol Neurosurg Psychiatry, 2011. 82(7): p. 761-5. Putzke, J.D., et al., Essential tremor: predictors of disease progression in a clinical cohort. J Neurol Neurosurg Psychiatry, 2006. 77(11): p. 1235-7. A. Schrag, et al., Essential tremor- an overdiagnosed condition? J Neurol 2000(247): p. 955-959. Samay Jain, Steven E. Lo, and E.D. Louis, Common Misdiagnosis of a Common Neurological Disorder. Arch Neurol., 2006. 63: p. 1100-1104. Elan D. Louis, et al., How Normal is ‘Normal’? Arch Neurol., 1998. 55: p. 222-227. Lawrenson, C., et al., The mystery of the cerebellum: clues from experimental and clinical observations. Cerebellum Ataxias, 2018. 5: p. 8. Filip, P., et al., Linking Essential Tremor to the Cerebellum: Physiological Evidence. Cerebellum, 2016. 15(6): p. 774-780. E.D. Louis, et al., Neuropathologic findings in essential tremor. Neurology, 2006. 13(66): p. 1756-9. Trujillo Diaz, D., et al., Banking brains: a pre-mortem "how to" guide to successful donation. Cell Tissue Bank, 2018. 19(4): p. 473-488. Choe, M., et al., Purkinje cell loss in essential tremor: Random sampling quantification and nearest neighbor analysis. Movement Disorders, 2016. 31(3): p. 393-401. Jordan E. Axelrad, et al., Reduced Purkinje Cell Number in Essential Tremor. Arch Neurol., 2008. 65(1): p. 101-107. Louis, E.D., et al., Torpedoes in Parkinson's disease, Alzheimer's disease, essential tremor, and control brains. Mov Disord, 2009. 24(11): p. 1600-5. Babij, R., et al., Purkinje cell axonal anatomy: quantifying morphometric changes in essential tremor versus control brains. Brain, 2013. 136(Pt 10): p. 3051-61. Yu, M., et al., Increased number of Purkinje cell dendritic swellings in essential tremor. Eur J Neurol, 2012. 19(4): p. 625-30. Louis, E.D., et al., Reduced Purkinje cell dendritic arborization and loss of dendritic spines in essential tremor. Brain, 2014. 137(Pt 12): p. 3142-8. Lee, D., et al., Climbing fiber-Purkinje cell synaptic pathology across essential tremor subtypes. Parkinsonism Relat Disord, 2018. 51: p. 24-29. Kuo, S.H., et al., Increased number of heterotopic Purkinje cells in essential tremor. J Neurol Neurosurg Psychiatry, 2011. 82(9): p. 1038-40.

ur na

1. 2.

17. 18. 19.

26. 27. 28. 29.

30.

31. 32.

33.

Jo

34.

of

25.

ro

24.

-p

23.

re

22.

lP

21.

Louis, E.D., et al., Heterotopic Purkinje Cells: a Comparative Postmortem Study of Essential Tremor and Spinocerebellar Ataxias 1, 2, 3, and 6. Cerebellum, 2018. 17(2): p. 104-110. Lee, P.J., et al., A Quantitative Study of Empty Baskets in Essential Tremor and Other Motor Neurodegenerative Diseases. J Neuropathol Exp Neurol, 2018. Cordelia R. Erickson-Davis, et al., ‘‘Hairy Baskets’’ Associated With Degenerative Purkinje Cell Changes in Essential Tremor. Neuropathol Exp Neurol, 2010. 69(3): p. 262-271. Lee1, M., et al., Decreased EAAT2 protein expression in the essential tremor cerebellar cortex. Acta Neuropathologica Communications, 2014. 2(157): p. 1-11. Louis E.D., et al., The Washington Heights-Inwood Genetic Study of Essential Tremor: Methodologic Issues in Essential-Tremor Research. Neuroepidemiology, 1997. 16(3): p. 124-133. Pan1, Y., et al., Whole tumor RNA-sequencing and deconvolution reveal a clinicallyprognostic PTEN:PI3K-regulated glioma transcriptional signature. Oncotarget, 2017. 8(32): p. 52474-52487. Love, M.I., W. Huber, and S. Anders, Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol, 2014. 15(12): p. 550. Pedregosa, F., et al., Scikit-learn- Machine Learning in Python. Journal of Machine Learning Research, 2011. 12: p. 2825-2830. Mathias Uhlén, et al., Tissue-based map of the human proteome. Science Proteomics, 2015. 347(6220): p. 394-405. Szklarczyk, D., et al., The STRING database in 2017: quality-controlled protein-protein association networks, made broadly accessible. Nucleic Acids Res, 2017. 45(D1): p. D362-D368. Wang, J., et al., WebGestalt 2017: a more comprehensive, powerful, flexible and interactive gene set enrichment analysis toolkit. Nucleic Acids Res, 2017. 45(W1): p. W130-W137. Fabregat, A., et al., The Reactome Pathway Knowledgebase. Nucleic Acids Res, 2018. 46(D1): p. D649-D655. Mi, H., et al., PANTHER version 11: expanded annotation data from Gene Ontology and Reactome pathways, and data analysis tool enhancements. Nucleic Acids Res, 2017. 45(D1): p. D183-D189. Mootha, V.K., et al., PGC-1α-responsive genes involved in oxidative phosphorylation are coordinately downregulated in human diabetes. Nature Genetics. 34(3): p. 267-273. Aravind Subramaniana, et al., Gene set enrichment analysis- A knowledge-based approach for interpreting genome-wide expression profiles. PNAS, 2005. 102(43): p. 15545-15550. Uhlén, M., et al., Tissue-based map of the human proteome. Science, 2015. 347(6220). Jimenez-Jimenez, F.J., et al., Update on genetics of essential tremor. Acta Neurol Scand, 2013. 128(6): p. 359-71. Louis, E.D., et al., Cerebellar Pathology in Familial vs. Sporadic Essential Tremor. The Cerebellum, 2017. 16(4): p. 786-791. Schmouth, J.F., P.A. Dion, and G.A. Rouleau, Genetics of essential tremor: from phenotype to genes, insights from both human and mouse studies. Prog Neurobiol, 2014. 119-120: p. 1-19.

ur na

20.

35. 36. 37. 38.

46. 47. 48. 49. 50. 51. 52. 53. 54. 55.

Jo

56.

of

45.

ro

43. 44.

-p

42.

re

41.

lP

40.

Muller, S.H., et al., Genome-wide association study in essential tremor identifies three new loci. Brain, 2016. 139(Pt 12): p. 3163-3169. Tio, M. and E.K. Tan, Genetics of essential tremor. Parkinsonism Relat Disord, 2016. 22 Suppl 1: p. S176-8. Deng, H., W. Le, and J. Jankovic, Genetics of essential tremor. Brain, 2007. 130(Pt 6): p. 1456-64. Liu, X., et al., Identification of candidate genes for familial early-onset essential tremor. Eur J Hum Genet, 2016. 24(7): p. 1009-15. Smith, P., et al., A Drosophila Model of Essential Tremor. Sci Rep, 2018. 8(1): p. 7664. Odgerel, Z., et al., Whole Genome Sequencing and Rare Variant Analysis in Essential Tremor Families. PLoS One, 2019. 14(8): p. e0220512. Louis, E.D., et al., Mapping Purkinje Cell Placement Along the Purkinje Cell Layer: an Analysis of Postmortem Tissue from Essential Tremor Patients vs. Controls. Cerebellum, 2016. 15(6): p. 726-731. Louis, E.D., et al., Structural study of Purkinje cell axonal torpedoes in essential tremor. Neurosci Lett, 2009. 450(3): p. 287-91. Baird, F.J. and C.L. Bennett, Microtubule defects & Neurodegeneration. J Genet Syndr Gene Ther, 2013. 4: p. 203. Gonatas, N.K., A. Stieber, and J.O. Gonatas, Fragmentation of the Golgi apparatus in neurodegenerative diseases and cell death. J Neurol Sci, 2006. 246(1-2): p. 21-30. Saya Nakagomia, et al., A Golgi fragmentation pathway in neurodegeneration. Neurobiol Dis, 2008. 29(2): p. 221-231. Alvarez-Miranda, E.A., M. Sinnl, and H. Farhan, Alteration of Golgi Structure by Stress: A Link to Neurodegeneration? Front Neurosci, 2015. 9: p. 435. Gruol, D., M. Manto, and D. Haines, Ca2+ signaling in cerebellar Purkinje neurons-editorial. Cerebellum, 2012. 11(3): p. 605-8. Hoxha, E., et al., Modulation, Plasticity and Pathophysiology of the Parallel FiberPurkinje Cell Synapse. Front Synaptic Neurosci, 2016. 8: p. 35. Handforth, A., Linking Essential Tremor to the Cerebellum-Animal Model Evidence. Cerebellum, 2016. 15(3): p. 285-98. Mark, M.D., et al., Keeping Our Calcium in Balance to Maintain Our Balance. Biochem Biophys Res Commun, 2017. 483(4): p. 1040-1050. Hisatsune, C., K. Hamada, and K. Mikoshiba, Ca(2+) signaling and spinocerebellar ataxia. Biochim Biophys Acta Mol Cell Res, 2018. 1865(11 Pt B): p. 1733-1744. Louis, E.D., 'Essential tremor' or 'the essential tremors': is this one disease or a family of diseases? Neuroepidemiology, 2014. 42(2): p. 81-9.

ur na

39.

Figure 1: Analysis of Variance in ET Cerebellar Transcripts A: Principal component analysis of differentially expressed gene transcripts for controls (blue) and ET (red) with two outliers removed. Principal components are the eigenvectors of the covariance matrix of the data, providing a broad overview of each sample’s gene expression profile. Clinical diagnosis was applied post PCA to explore ET versus control differences. B: Mean PCA values for ET and control samples with SEM. Two tailed T-Test comparing

ur na

lP

re

-p

ro

of

Control to ET for PC2 * p=0.020.

Figure 2: String Network Clusters of Dysregulated Genes in ET

Jo

A. Results from String Network database for 230/231 ‘top gene hits’. One transcript did not map, LINC00641 (non-coding RNA gene). Clustered into the top 7 (kmean) groups show significant clustering with axon guidance (Red), kinesin/microtubule (Dark Blue), ribosomal (Yellow), p53/cell cycle (Light Blue), Golgi-ER transport (Green), calcium signaling (Purple); all others are in Tan. B. String Networks for highly relevant pathways for axon guidance (Red), kinesins/microtubule motor activity (Dark Blue), ER-Golgi transport (Green), and calcium signaling (Purple).

of ro -p re

Figure 3: ET Clusters of Pathway Dysregulation

lP

A. Heatmap from single sample GSEA (ssGSEA) for ET samples showing gene set enrichment for the dysregulated GO biological processes listed in Table 2. ET samples were hierarchically clustered based on normalized enrichment scores from the ssGSEA and grouped based on the

ur na

dendrogram. Color scale of red to green indicates high levels of enrichment to low levels of enrichment, respectively. Kruskal-Wallis H test p=0.0001. B. Spearman correlation matrix of all ET samples generated in Prism from the ssGSEA NES. Color scale of blue to red indicate spearman (r) correlation coefficients from positive to negative,

Jo

respectively. White indicates no correlation.

of

ro

-p

re

lP

ur na

Jo

Table 1: Summary of Dysregulated Processes from 231 ’Top Gene Hits’ Database Mining A compilation of dysregulated pathways identified via multiple curated databases based on differential gene expression analysis of ET versus control cerebellar samples. See methods/citations for information regarding each database. Number of hits indicates how many genes from the 231 top hits are involved in each pathway. False Discovery Rate (FDR) was

of

calculated via the database, with input of statistics from DEseq2.

Hits

Genes

Axon development

19

MEGF8, BMPR2, FZD3, DYNC1H1, KIF5A, NREP, ZIC1, BTBD3, DCAF8, TTBK2, SLITRK4, ENAH, FYN, RAPGEF2, MYCBP2, S100B, EPHA4, KIF2A, CAMSAP2

1.79e-13

Axonogenesis

12

EPHA4, SSH1, BMPR2, APOE, MEGF8, ADNP, ROCK2, CREB1, MAPK1, DMN3, SPG7, SETX

2.01e-10

Axon guidance

14

MEGF8, BMPR2, EPHA4, ARHGEF12, FYN, SSH1, ENAH, MAPK1, FZD3, PIK3R3, ROCK2, SLITRK4, ADNP, S100B

Cellular component organization

30

FDR

2.85e-08

lP

re

-p

ro

Biological Process

0.000147

Post-Translational Protein Modification

13

DYNC1H1, GFPT1, HNRNPC, ARCN1, MIA3, OPCML, ST8SIA3, WAC, RAB6B, RANBP2, USP7, CALM1, COPB2

0.00721

Synaptic transmission

9

ADCY1, ADNP, PRKCE, SYNJ1, RAB3GAP1, S100B, KIF5B, MAPK1, ATP2B1

0.00482

Regulation of synapse structure/activity

11

AGT, CREB1, MAPK1, APOE, YWHAG, SLITRK4, RAPGEF2, EEF2K, RAB3GAP1, ADNP, S100B

0.00657

Calcium Signaling

22

RYR1, FZD3, PRKCE, ADCY1, FAT2, RAB11F1P4, PRKCB, UNC80, CALD1, EEF2K, CALN1, MEGF8, RAPGEF2, ATP2B1, CAPN2, MYO6, CERK, AHCYL1, CALM1, S100B, CDS2, SPTAN1

0.004

Modulation of synaptic transmission

14

CPLX3, AGT, CREB1, NISCH, MAPK1, PRKCE, APOE, YWHAG, RAPGEF2, SYNJ1, RAB3GAP1, ADNP, S100B, KIF5B

0.0135

Jo

ur na

AGT, CAND1, JUN, BMPR2, CREB, MAPK1, PRKCE, APOE, DNM3, YWHAG, RPS3, KIAA0947, JARID2, SLITRK4, ROCK2, MEGF8, FYN, RAPGEF2, CD47, SYNJ1, EEF2K, RAB3GAP1, ADNP, NCKAP1, HUWE1, FNIP1, EPHA4, TSC1, SETX, KIF5B

KIF2A, KIF5B, KIF5A, COPB2, DYNC1H1, KLC2, KLC1, RAB6B, RAB3GAP1, COPB2

0.0197

Hits

Genes

FDR

Kinesin complex

5

KIF2A, KIF5A, KIF5B, KLC1, KLC2

0.000447

Microtubule associated complex

8

DYNC1H1, KIF5A, RANBP2, KLC1, SBTBN5, KLC2, KIF2A, KIF5B

0.000955

Golgi apparatus

12

GOLGA8A, TGOLN2, FUT9, ST8SIA3, RAB6B, AFTPH, COPB2, ARCN1, NUMA1, SCAMP1, NMNAT2, CPD

0.0149

Golgi membrane

8

FUT9, TGOLN2, AFTPH, ST8SIA3, RAB6B, CALN1, NUMA1, MGAT5

0.0162

Vesicle

20

DYNC1H1, KIF5A, TGOLN2, KLC1, RAB6B, MIA3, EXOC5, COPB2, ARNC1, SPTBN5, KLC2, KIF2A, SCAMP1, SPTAN1, CPLX3, DNM3, MYO6, PRKCE, APOE, CD47

0.0215

Molecular Function

ro

Cellular Component

of

10

-p

Retrograde vesiclemediated transport, Golgi-ER

KEGG Pathways

6

DYNC1H1, KIF5A, KLC1, KLC2, KIF2A, KIF5B

0.000211

Hits

Genes

FDR

FZD3, PIK3R3, BMPR2, EPHA4, SSH1, MAPK1, ROCK2, ARHGEF12, ENAH, FYN

0.000674

10

Jo

ur na

Axon Guidance

FDR

lP

Microtubule motor activity

Genes

re

Hits

of

Table 2: GSEA Ontology Pathways

Representative enriched pathways from GSEA of PCA values. Pathways correspond to correlated (positive) and anti-correlated

ro

(negative) PCA space with respect to Biological Processes, Cellular Component and Molecular function. Enrichment scores calculated via GSEA software. To calculate enrichment score, GSEA begins with the top of the ranked gene list and adds to the running sum if a

-p

gene is a member of the candidate gene set. If it is not, the number subtracts. This process is repeated for each gene in the ranked list.

ur na

GO_PROTEIN_LOCALIZATION_TO_ENDOPLASMI C_RETICULUM GO_ESTABLISHMENT_OF_PROTEIN_LOCALIZATI ON_TO_ENDOPLASMIC_RETICULUM GO_TRANSLATIONAL_INITIATION GO_ESTABLISHMENT_OF_PROTEIN_LOCALIZATI ON_TO_ORGANELLE GO_PROTEIN_TARGETING_TO_MEMBRANE

Jo

GO_PROTEIN_TARGETING GO_ESTABLISHMENT_OF_PROTEIN_LOCALIZATI ON_TO_MEMBRANE GO_PROTEIN_LOCALIZATION_TO_MEMBRANE GO_TRANSLATIONAL_ELONGATION GO_TRANSLATIONAL_TERMINATION

Cellular Component

lP

Biological Process

Sc or e

re

The enrichment score for that gene set is equal to the largest absolute value that the running sum achieves.

7.7 59 7.7 00 7.6 44 6.6 00 6.3 79 6.0 94 5.7 49 5.5 47 4.7 27 4.6

GO_TRANSPORT_VESICLE GO_MYELIN_SHEATH GO_PRESYNAPSE GO_VESICLE_MEMBRANE GO_DENDRITE GO_TRANSPORT_VESICLE_M EMBRANE GO_POSTSYNAPSE GO_AXON GO_EXOCYTIC_VESICLE GO_CELL_BODY

Sc or e 5.8 21 5.7 24 5.5 60 5.1 02 4.7 20 4.7 04 4.6 74 4.6 40 4.5 96 4.3

Molecular Function GO_CLATHRIN_BINDING GO_COENZYME_BINDING GO_COFACTOR_BINDING GO_GDP_BINDING GO_GUANYL_NUCLEOTIDE_BINDIN G GO_NAD_BINDING GO_TRANSLATION_FACTOR_ACTIVI TY_RNA_BINDING GO_UBIQUITIN_LIKE_PROTEIN_LIG ASE_BINDING GO_GTPASE_ACTIVITY GO_TRANSLATION_INITIATION_FAC

Sc or e 3.1 76 2.8 94 3.4 63 3.6 08 4.0 88 2.9 53 3.7 81 3.0 77 4.9 51 3.7

ur na

GO_ORGANELLE_LOCALIZATION

GO_CYTOPLASMIC_TRANSLATION GO_SYNAPTIC_SIGNALING

Jo

GO_VESICLE_TARGETING GO_ENDOPLASMIC_RETICULUM_ORGANIZATIO N

GO_EXOCYTOSIS GO_TRANSCRIPTION_COUPLED_NUCLEOTIDE_E XCISION_REPAIR GO_REGULATION_OF_VESICLE_MEDIATED_TRA

TOR_ACTIVITY GO_PROTEIN_TRANSPORTER_ACTI VITY GO_UBIQUITIN_LIKE_PROTEIN_TRA NSFERASE_ACTIVITY

of

GO_SYNAPTIC_MEMBRANE GO_NEURONAL_POSTSYNAP TIC_DENSITY

GO_TERMINAL_BOUTON GO_ER_TO_GOLGI_TRANSPO RT_VESICLE_MEMBRANE GO_POSTSYNAPTIC_MEMBR ANE

lP

GO_REGULATION_OF_SYNAPTIC_PLASTICITY GO_REGULATION_OF_RESPONSE_TO_EXTRACE LLULAR_STIMULUS GO_REGULATION_OF_SYNAPSE_STRUCTURE_O R_ACTIVITY

51 3.8 28 3.6 49 3.4 93 3.4 61 3.3 57 3.1 48 2.9 93 2.9 34 2.9 02 2.8 99 2.8 97 2.8 84 2.8 79 2.8 68 2.8 50

ro

GO_VESICLE_LOCALIZATION GO_ER_TO_GOLGI_VESICLE_MEDIATED_TRANS PORT GO_POSTTRANSCRIPTIONAL_REGULATION_OF_ GENE_EXPRESSION

GO_EXCITATORY_SYNAPSE

-p

GO_VESICLE_ORGANIZATION

GO_AXON_PART GO_GOLGI_ASSOCIATED_VES ICLE

re

GO_GOLGI_VESICLE_TRANSPORT GO_MODULATION_OF_SYNAPTIC_TRANSMISSI ON

44 4.5 96 4.5 26 4.5 01 4.4 97 4.2 91 4.2 89 4.0 31 4.0 25 3.8 74 3.8 25 3.6 97 3.5 53 3.5 47 3.3 47 3.3 22 3.3 11 3.2

GO_SECRETORY_VESICLE GO_PRESYNAPTIC_MEMBRA NE GO_SITE_OF_POLARIZED_GR OWTH GO_MICROBODY_MEMBRAN E GO_PERIKARYON GO_ER_TO_GOLGI_TRANSPO RT_VESICLE GO_CLATHRIN_COAT_OF_CO ATED_PIT

26 3.5 73 3.2 28

of ro

GO_NEURAL_NUCLEUS_DEVELOPMENT GO_REGULATION_OF_NEUROTRANSMITTER_LE VELS

97 3.2 76 3.1 79 3.1 79

-p

NSPORT GO_REGULATION_OF_SYNAPTIC_VESICLE_TRA NSPORT

Sc or e

Biological Process

3.4 64

Cellular Component

Jo

ur na

lP

GO_SENSORY_PERCEPTION_OF_CHEMICAL_STI MULUS

re

PC1 Anti-Correlated

Biological Process

GO_PROTEIN_FOLDING

GO_RESPONSE_TO_GROWTH_FACTOR GO_NEGATIVE_REGULATION_OF_PHOSPHORYL

GO_CENTROSOME

GO_CENTRIOLE GO_MICROTUBULE_ORGANI ZING_CENTER_PART

Sc or e 3.1 81 3.1 65 3.0 95

Molecular Function

Sc or e

GO_RHO_GUANYL_NUCLEOTIDE_E XCHANGE_FACTOR_ACTIVITY

3.3 54 3.1 71 3.0 50

Molecular Function

Sc or e

GO_PROTEIN_SERINE_THREONINE_ KINASE_ACTIVITY GO_OLFACTORY_RECEPTOR_ACTIVI TY

PC2 Correlated Sc or e 4.5 28 4.3 61 4.2

Cellular Component GO_AXON GO_MYELIN_SHEATH GO_DENDRITE

Sc or e 3.8 55 3.8 49 3.7

GO_UBIQUITIN_LIKE_PROTEIN_LIG ASE_BINDING GO_GTPASE_ACTIVITY GO_UNFOLDED_PROTEIN_BINDING

4.4 63 3.6 96 3.6

ro

GO_CELL_BODY

ur na

GO_PROTEIN_STABILIZATION

Jo

GO_NEURON_PROJECTION_MORPHOGENESIS GO_NEGATIVE_REGULATION_OF_INTRACELLUL AR_SIGNAL_TRANSDUCTION GO_REGULATION_OF_TRANSMEMBRANE_TRA NSPORT GO_TRANSMEMBRANE_RECEPTOR_PROTEIN_T YROSINE_KINASE_SIGNALING_PATHWAY GO_NEGATIVE_REGULATION_OF_PROTEIN_SER INE_THREONINE_KINASE_ACTIVITY GO_SIGNAL_TRANSDUCTION_BY_PROTEIN_PH

04 3.6 99 3.6 84 3.5 96 3.3 50 3.2 58 3.2 40 3.2 13 3.1 08 3.0 33 3.0 07 2.9 99 2.9 65 2.9 09 2.9 04 2.9 00 2.8 93

GO_POSTSYNAPSE

-p

GO_ANCHORING_JUNCTION GO_PRESYNAPSE GO_AXON_PART

GO_EXCITATORY_SYNAPSE GO_CELL_SUBSTRATE_JUNCT ION GO_MEMBRANE_MICRODO MAIN GO_TRANSPORT_VESICLE GO_NUCLEAR_BODY GO_CELL_CELL_ADHERENS_J UNCTION GO_EXOCYTIC_VESICLE_ME MBRANE GO_UBIQUITIN_LIGASE_COM PLEX GO_SYNAPTIC_MEMBRANE GO_MICROTUBULE

GO_GUANYL_NUCLEOTIDE_BINDIN G GO_TRANSCRIPTION_FACTOR_BIND ING GO_UBIQUITIN_LIKE_PROTEIN_TRA NSFERASE_ACTIVITY

of

GO_PHOSPHATASE_COMPLE X

lP

GO_SYNAPTIC_SIGNALING GO_RESPONSE_TO_TOPOLOGICALLY_INCORRE CT_PROTEIN GO_NEGATIVE_REGULATION_OF_KINASE_ACTI VITY GO_NEGATIVE_REGULATION_OF_TRANSFERAS E_ACTIVITY GO_MODULATION_OF_SYNAPTIC_TRANSMISSI ON GO_RESPONSE_TO_ENDOPLASMIC_RETICULUM _STRESS GO_CELL_MORPHOGENESIS_INVOLVED_IN_NE URON_DIFFERENTIATION GO_SMALL_GTPASE_MEDIATED_SIGNAL_TRAN SDUCTION GO_REGULATION_OF_INTRINSIC_APOPTOTIC_S IGNALING_PATHWAY

39 4.1 10 4.0 32 3.9 49 3.9 05 3.8 67 3.7 88 3.6 66 3.5 93 3.5 41 3.5 16 3.4 52 3.4 42 3.4 08 3.3 87 3.3 12 3.3 12 3.3

re

ATION GO_REGULATION_OF_APOPTOTIC_SIGNALING_ PATHWAY

GO_CHROMATIN_BINDING GO_ATPASE_REGULATOR_ACTIVITY GO_GDP_BINDING GO_HEAT_SHOCK_PROTEIN_BINDI NG GO_G_PROTEIN_COUPLED_RECEPT OR_BINDING GO_UBIQUITIN_LIKE_PROTEIN_CO NJUGATING_ENZYME_ACTIVITY GO_TRANSCRIPTION_COACTIVATO R_ACTIVITY

69 3.6 21 3.6 00 3.4 50 3.3 30 3.1 99 3.1 95 3.1 67 3.0 11 2.9 80 2.9 69

GO_NEURON_PROJECTION_GUIDANCE

ur na

GO_REGULATION_OF_CELL_GROWTH

GO_PROTEIN_LOCALIZATION_TO_NUCLEUS GO_REGULATION_OF_MICROTUBULE_BASED_P ROCESS

Jo

GO_REGULATION_OF_EXOCYTOSIS

Biological Process

GO_CILIUM_ORGANIZATION

of ro

lP

GO_REGULATION_OF_PROTEIN_STABILITY GO_CELLULAR_RESPONSE_TO_TOPOLOGICALLY _INCORRECT_PROTEIN

-p

GO_NUCLEAR_TRANSPORT GO_NEGATIVE_REGULATION_OF_APOPTOTIC_S IGNALING_PATHWAY GO_POSITIVE_REGULATION_OF_CELL_PROJECTI ON_ORGANIZATION

05 3.3 03 3.2 48 3.2 34 3.2 26 3.2 01 3.1 02 3.0 80 3.0 64 3.0 30 2.9 97 2.9 84 2.9 73

re

OSPHORYLATION GO_REGULATION_OF_SYNAPSE_STRUCTURE_O R_ACTIVITY GO_POSITIVE_REGULATION_OF_NEURON_DIFF ERENTIATION

PC2 Anti-Correlated Sc or e 3.4 12

Cellular Component

Sc or e

Molecular Function

Sc or e

GO_INTRACILIARY_TRANSPO RT_PARTICLE

2.9 10

GO_TRNA_METHYLTRANSFERASE_A CTIVITY

3.2 96

ro

GO_CILIUM_MORPHOGENESIS

of

3.3 45

-p

Supplemental Table 1: Demographics and Postmortem Features with RNA Quality Measures Table includes sample identifiers, age at death, sex, diagnosis, brain weight in grams, Braak staging, CERAD staging, post-mortem

re

interval to freezing (PMI-Frozen), RNA integrity number (RIN), percent alignment to the genome and the brain bank location it was received from.

lP

CT = control, ET = essential tremor, CERAD = Consortium to Establish a Registry for Alzheimer’s disease, RIN = RNA Integrity

ur na

Number (maximum = 10), ETCBR = Essential Tremor Centralized Brain Repository, NA = Not Available

Supplemental Table 2: 231 ‘Top Gene Hits’ in ET Versus Control Cerebellar Cortex Transcriptome List of the top 231 differentially expressed gene transcripts with a false discovery rate (FDR) of less than 0.25 (25%). Differential

Jo

gene expression data from DESeq2. A positive log2FoldChange indicates increased expression in ET versus controls, and vice versa.

Supplemental Table 3: Heterogeneity of ET Cerebellar Cortex Transcripts Custom created heatmap evaluating the current knowledge available on 231 ‘top gene hits’ in the field of Essential Tremor. Can sort heatmap by any column, all columns contain values corresponding to the color code.

Columns:

of

1. Gene name 2. Log2 Fold Change

ro

3. pValue

5. PC1 Correlations 6. PC2 Correlations

re

7. Purkinje cell (PC) staining (Human protein atlas)

-p

4. pAdjusted (FDR)

8. Molecular layer (ML) staining (Human protein atlas)

lP

9. Granule cell layer (GCL) staining (Human protein atlas)

10. Reads per kilobase of transcript, per million (RPKM) Cerebellar Hemisphere (Human protein atlas) 11. RPKM Cerebellum (Human protein atlas)

ur na

12. PubMed: Cerebellum + gene name 13. PubMed: GABA + gene name

14. PubMed: Purkinje Cell + gene name 15. PubMed: Brain + gene name

Jo

16. PubMed: Total publication for gene

Supplemental Figure 1: Outlier Identification via PCA Principal component analysis (PCA) of differentially expressed gene transcripts for controls (blue) and ET (red). Two outlier samples (circled in black) were found to be under-sequenced and removed from analysis.

Supplemental Figure 2: Pathological Correlates with PCA

of

Correlation of LH&E Torpedo counts (ovoid axonal swellings) per Purkinje cell (PC) with principal component 2 (PC2). (Spearman

Jo

ur na

lP

re

-p

ro

correlation p=0.0049). Control (blue), ET (red).