Single nucleus sequencing of mouse hippocampi four days after pilocarpine induced status epilepticus

All current drug treatments for epilepsy, a neurological disorder affecting over 50 million people merely treat symptoms, and a third of patients with epilepsy do not respond to medication. There are no disease modifying treatments that may be administered briefly to patients to enduringly eliminate spontaneous seizures and reverse cognitive deficits. Applying network approaches to whole tissue and single-nuclei transcriptomic data collected from mouse models of temporal lobe epilepsy and publicly available transcriptomic data from human temporal lobectomy samples, we confirmed a previously described pattern of rapid and transient induction of the Janus kinase/signal transducers and activators of transcription (JAK/STAT) pathway within days of epileptogenic insult. This was followed by a resurgent activation weeks to months later with the onset of spontaneous seizures. Targeting the first wave of JAK/STAT activation after epileptic insult did not prevent seizures. However, inhibition of the second wave with CP690550 (Tofacitinib) over a 2-week period enduringly suppressed seizures, rescued deficits in spatial memory, and alleviated epilepsy-associated histopathological alterations. Seizure suppression lasted for at least 2 months after the final dose. These results indicate that reignition of inflammatory JAK/STAT3 signaling in chronic epilepsy opens a window for disease modification with the FDA-approved, orally available drug CP690550.

Methods
Pilocarpine model of status epilepticus. Status epilepticus was induced in male and female C57BL/6CRN mice with pilocarpine (8-12 weeks old). Mice were injected with terbutaline (4 mg/kg i.p., Sigma T2528) and methylscopolamine (4 mg/kg ip, Sigma S8502) to alleviate respiratory and cardiovascular effects of pilocarpine, then 30 min later received saline or pilocarpine (280 mg/kg i.p. as free base, Sigma P6503). After 1 hour of SE, seizures were interrupted by diazepam (10mg/kg i.p.). Mice were weighed and scored for behavioral deficits (Irwin score) daily, with injections of lactated Ringer’s solution as needed.

Tissue Harvest and Sequencing. To select which epileptic mice would be sent for sequencing, we determined a composite score of seizure severity during pilocarpine SE, weight loss 24 hours after SE, ability to build a nest and neurobehavioral recovery scores on day 4 after SE (See Fig. S14 in linked manuscript Hoffman et al). The mean number of standard deviations from the mean of each of these 4 measures was calculated across all 6 saline-treated and 14 pilocarpine-treated mice, and mice were selected based on wide separation between saline and pilo groups. Both hippocampi and a small portion of neocortex overlaying one hippocampus were snap frozen and shipped to Novagene.

Sequencing, Preprocessing, Expression Quantification and Quality Control. Each sample was individually subjected to library preparation and Illumina sequencing. Reads were processed and mapped with Cell Ranger software (10x Genomics) and the Seurat R package at Novagene.  Counts matrices for each sample were combined into an adata object in the Roopra lab and subjected to quality control, cell and gene filtering using the Scanpy package. Genes were kept that were expressed in at least 100 cells across all samples. Cells with genes numbering between 1000 and 6500 and with fewer than 5% mitochondrial genes were kept. This resulted in 64,647 cells and 21,423 genes.

For PCA analysis of samples, reads were summed across all cells for each sample to produce a pseudo-bulked sample (n=10) x gene (21,423) array.  4000 most variant genes were selected to subject samples to PCA analysis.  All 5 pilo samples clustered together when projecting either the 1st 2 or the 2nd 2  PCA components. Four saline samples clustered together away from the pilo samples.  One saline sample did not cluster with either group and was discarded (Fig. S14 in linked manuscript Hoffman et al).  This resulted in 5 pilo samples, 4 saline samples, and 58,365 cells.

Normalization and Clustering. Data was normalized using the scanpy preprocessing normalization tool.  Highly expressed genes were omitted, and the target sum was set to 10,000 (transcripts per 10k or TPTK). Counts were converted to log(1+tptk) and highly variable genes identified with min_mean=0.0125, max_mean=3, min_disp=0.5. Genes were scaled to unit variance and values>10 standard deviations were clipped. PCA with 50 components was performed and used to generate a neighborhood graph with n_neighbors=15. Leiden clustering with resolution=0.25 yielded 27 clusters.

Annotations. 

Class annotation: Clusters were allocated to classes using the following markers: Inhibitory Neurons: GAD1, Excitatory Neurons: SLC17A7, Astrocytes: AQP4, Oligodendrocytes: SOX10, Endothelial cells: RGS5, Myeloid: PTPRC.  Each Class were broken up into Cell Classes as follows:

Excitatory Neurons: The single cell dataset from Yao et al. (ref 5 in linked manuscript Hoffman et al) was converted to an anndata object and trimmed to glutamatergic neurons (CA1-ProS, CA2-IG-FC, CA3, CR, DG). DEGs/markers were identified for each subtype using the rank_genes_groups scanpy function with method=’wilcoxon’. The same was done for our Excitatory adata object using groupby=’leiden’.  The top 100 Yao markers per cell type and our leiden clusters were used to perform pairwise fisher exact tests to look for overlaps between Yao cell types and our leiden clusters. The Yao cell type with highest overlap with Leiden pairs as judged by highest score (log(odds)*-log(p)) was assigned to the cluster. This yielded the following subclasses: CA1_ProS, CA1_ve, CA1_do, CA3, CR, DG. Using Hipposeq (6 in linked manuscript Hoffman et al) the CA3 subclass (markers: [SPOCK1, MGAT4C, ELAVL2]) was broken up into CA3_ve, CA3_do and mossy cells using the following markers: CA3_ve: [NECAB1,COCH,ADGRA1,CPNE7], CA3_do: [CHGB,RIMBP2], mossy: [CALCRL,CALB2]. Yao et al CR cell markers were RELN, CACNA2D2, STMN1. DG markers were PROX1, PCP4, ADARB2, STXBP6.

Inhibitory Neurons: As with excitatory neurons, the Yao dataset was trimmed to GABAergic neurons. Rank genes groups was run on the inhibitory adata object and the Yao adata object. Fisher exact tests with the top 100 genes per group was performed for each Yao subtype and Inhibitory leiden cluster. This approach defined 10 inhibitory subtypes in our dataset: Lamp5 Lhx6, Lamp5, Meis2, Ntng1, Pvalb Vipr2, Pvalb, Sst Ntng1, Sst, Vip, Vip Igfbp6.

Astrocytes: AQP4 was used to define astrocytes. GFAP was used to define activated and resting states.

Myeloid cells: The myeloid cell class was first defined by PTPRC. The myeloid adata object was then reclustered with n_neighbors=15 and leiden resolution=0.25. This gave 12 clusters. GFAP defined activated microglia. SIGLCECH positive and PTPRC negative marked resting microglia. Cells with high TOP2A, PCNA and MKI67 defined proliferating microglia. CCR2 and CD44 positive denoted invading monocytes. MRC1 and VCAM marked Perivascular Macrophages.

Oligodendrocytes: The oligodendrocyte cluster was reclustered with n_neighbors=15 and leiden resolution=0.25 to give 9 clusters. Mature oligodendrocytes were marked by MOG. Progenitors (OPC in Fig. 2 in linked manuscript Hoffman et al) were marked by VCAN.