A 4D transcriptomic map for the evolution of multiple sclerosis-like lesions in the marmoset brain

Single-time-point histopathological studies on postmortem multiple sclerosis (MS) tissue fail to capture lesion evolution dynamics, posing challenges for therapy development targeting development and repair of focal inflammatory demyelination. To close this gap, we studied experimental autoimmune encephalitis (EAE) in the common marmoset, the most faithful animal model of these processes. Using MRI-informed RNA profiling, we analyzed ~600,000 single-nucleus and ~55,000 spatial transcriptomes, comparing them against EAE inoculation status, longitudinal radiological signals, and histopathological features. We categorized 5 groups of microenvironments pertinent to neural function, immune and glial responses, tissue destruction and repair, and regulatory network at brain borders. Exploring perilesional microenvironment diversity, we uncovered central roles of EAE-associated astrocytes, oligodendrocyte precursor cells, and ependyma in lesion formation and resolution. We pinpointed imaging and molecular features capturing the pathological trajectory of WM, offering potential for assessing treatment outcomes using marmoset as a platform.


Sample index
105 and is also made available for use under a CC0 license.
(which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25.559371doi: bioRxiv preprint FigS1.related to Fig1.Experimental design used to create a transcriptome map for the evolution of white matter (WM) lesions with single nucleus resolution.
(A) The dataset includes marmosets that were inoculated with human white matter (hWM) or recombinant human myelin oligodendrocyte glycoprotein (hMOG) emulsified in complete (CFA) or incomplete Freund's adjuvant (IFA) to induce experimental autoimmune encephalomyelitis (EAE) and healthy, naïve controls.
(B) The experimental workflow involved scanning and categorizing brain tissue using MRI.Postmortem brains were sliced into 3-mm slabs from anterior (A) to posterior (P).Specific areas of interest were sampled as cylinders with a diameter of 2 mm and height of 3 mm.These sampled areas were labeled on the standard slab (SS) index and grouped based on disease condition.The number of biological repeats (BR) is indicated by black annotations on the SS.
(C) Nuclei were isolated from the sampled areas to prepare cDNA libraries, which were then sequenced.
105 and is also made available for use under a CC0 license.
(which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC  (A) The experimental workflow involves identifying WM lesions using magnetic resonance imaging (MRI) and preparing postmortem tissue for spatial transcriptome analysis using the 10x Visium platform.
(C) Proton density-weighted (PDw) MRI was performed at the terminal time point and matched with the sampled tissue areas.
(D) A postmortem MRI atlas overlaid with brain region labels was used to match with the selected terminal PDw MRI.
105 and is also made available for use under a CC0 license.
(which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25.559371doi: bioRxiv preprint  (which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25.559371doi: bioRxiv preprint FigS3.related to Fig1.Matched brain regions imaged by different modalities.i) In vivo brain PDw MRI acquired on a 7 Tesla scanner at the disease terminal.
ii) Histological examination of myelin content using Sudan black (SB) and nuclear fast red (NFR) staining of postmortem tissue.
iii) RNA landscape visualized through the 10x Visium platform.iv) Marmoset MRI atlas with region annotations.
See source data for the full list of abbreviations for brain regions.
FigS4.related to Fig1.Image processing and resolution enhancement pipelines used to quantify white matter (WM) lesion load, generate lesion subregion masks, and deconvolute mixed transcript signals.
(A) To gain insight into lesion development and progression, an MRI characterization of MS-like lesions in the Marmoset Quantitatively (M3Q) pipeline was developed, including the following steps (Methods): Proton density-weighted (PDw) MRI images at baseline (before EAE induction) and terminal (before tissue collection) time points were subjected to the N4 bias field correction algorithm.The brain portion of the images, corresponding to the region of interest, was extracted using a skull-removal algorithm to improve image alignment.Each image was individually registered to the marmoset MRI atlas using bUnwarpJ, achieving spatial alignment across time points and animals.MRI intensity changes were calculated by subtracting the normalized terminal image from the baseline image, providing information about lesion location.Binary lesion masks were created by applying intensity thresholding to the subtracted images, segmenting the region of interest.The WM portion of the lesion mask was parsed and analyzed using atlas annotation indexing, enabling further characterization of the lesions based on location and distribution.
(B) To gain insight into the regionally enriched signal distribution within and near the WM lesion, a spatial transcriptome (ST) image processing pipeline was employed with the following steps: The myelinated WM area (Sudan black-positive) was extracted using the "Color Deconvolution" function in Fiji with the default "H DAB" setting, resulting in an SB + WM binary mask.The coordinates of the SB + WM mask were transferred to the 10x Visium spot hexagon coordinate system using Seurat, facilitating the spatial mapping of gene expression data.To distinguish SB -gray matter (GM) from SB -demyelinated WM, GM and lesion gene module scores were calculated and filtered to create GM and lesion masks accordingly.Spots that exhibited both SB -and IMM + signals were identified as the lesion core.Next, 10 concentric rims (SB + WM_rims) extending outward from the lesion core were assigned to mark the adjacent lesion neighborhoods.The normal-appearing (NA) WM area was annotated by subtracting the lesion neighborhoods from the SB + WM mask in animals with experimental autoimmune encephalomyelitis (EAE), and this region was labeled as "SB + WM_NA.Ctrl."Additionally, lesion neighborhoods that overlapped with the GM mask were labeled as "SB -notWM_rims," while the supplemental area was labeled as "SB -notWM_EAE" in animals with EAE.
Subregions within the lesion core were further divided based on centripetal rim assignments (SB -WM_rims).For healthy animals, "SB + WM_He.Ctrl" and "SB -notWM_He" labels were used to annotate tissue with or without SB staining, respectively.
(C) The spatial transcriptome resolution, initially measured at the spot level using the 10x Visium platform, was further enhanced to the subspot level using the BayesSpace algorithm.The enhanced signals at the subspot level were then rescaled to a ratio of 1 across genes.To gain insights into cell-type distributions, the averaged expression of gene sets enriched in specific cell types, acquired from single-nucleus RNA sequencing (snRNA-seq) references, was calculated.Cell-type locations were inferred by assessing the relative profile similarity score, allowing identification and mapping of different cell types within the spatial transcriptome data.By employing these techniques, the spatial transcriptome analysis achieved higher resolution, enabling a more detailed understanding of the cellular composition and gene expression patterns within the tissue of interest.(A) To demonstrate the phenotypic characterization of the ME within matched brain regions, contrasts between: i) myelin content visualized through Sudan black (SB) and nuclear fast red (NFR) staining, ii) unbiased ME clustering achieved through transcriptome similarity analysis, iii) 5 ME groups assigned by ME profile similarity, and iv) lesion subregions assigned using rim analysis and overlaid onto SB/NFR stained images, were indexed.The "Color Deconvolution" (FigS4B) for Samples 1-4 was unsuccessful due to suboptimal contrast between SB and NFR staining, resulting in their exclusion from the lesion subregion assignment in the rim analysis; however, they are included for ME clustering analysis.Significantly (FDR < 0.05 & abs(Log2FC) > 0.5) enriched ME between pairs of subregional WM area are colored accordingly."IL.WM" contains SB-WM_-rim5 to SB-WM_-rim2, "PL.WM" contains SB-WM_-rim1 to SB+WM_rim1, "EL.WM" contains SB+WM_rim2 to SB+WM_rim10), "NA.WM" contains SB+WM_NA.Ctrl, "He.WM" contains SB+WM_He.Ctrl, "He.notWM" contains SB-notWM_He, "NA.notWM" contains SB-notWM_EAE, and "EL.noWM" contains SB-notWM_rims.
105 and is also made available for use under a CC0 license.
(which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC                                                 (A) UMAP scatter plots color-coded by microenvironment (ME) clustering and subregional labeling.The number of spots per microenvironment is indicated in parentheses.Abbreviations: SB (Sudan black), WM (white matter), NA (normal appearing), He (healthy), Ctrl (control).
Gene names starting with "*" indicate human (hs) or mouse (mm) orthologs of marmoset gene identification numbers (See Table S9 for the full list).
(B) Dot plot showing the averaged and scaled expression of selected genes across L2 subclusters.Genes are split into groups to aid label tracking.Abbreviations: ME marker (genes used to annotate ME groups in Fig2A), ME DEG (differentially expressed genes across 28 ME), rDEG (regional differentially expressed genes), ferro.(ferroptosis genes), compl.(complement genes), dOPC (differentiating OPC enriched genes).105 and is also made available for use under a CC0 license.
(which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25.559371doi: bioRxiv preprint    (A) -(G) The level 1 (L1) analysis identified 6 major cell classes, which were further divided into 7 partitions in the level 2 (L2) analysis.UMAP scatter plots show subclustering of the following cell classes: microglia (MIC), oligodendrocyte progenitor cells (OPC), astrocytes (AST), oligodendrocytes (OLI), peripherally derived immune cells (P.IMM), neurons (NEU), and vascular/meningeal/ventricular cells (VAS).The number of nuclei analyzed in each L2 UMAP plot is listed in parentheses.
105 and is also made available for use under a CC0 license.
(which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC      (A) The inner pie charts show the relative nuclei proportion of L1 cell classes, and the outer donut charts display L2 sub-clusters.Both control and EAE samples for each tissue type are included, with the total number of nuclei listed in parentheses.In control animals, the composition of L1 cell classes varied across tissue types, as expected.Specifically, a higher number of glial cells were found in parietal white matter (pWM), while neurons were predominant in parietal cortex (pCTX) and lateral geniculate nucleus (LGN) region.In EAE animals, there was a significant expansion of microglia (MIC) and peripheral immune cells (P.IMM) partitions in all tissue types.
(B) Donut charts show the relative nuclei proportion of glial and immune clusters in EAE animals across different tissue types.The compositions of the P.IMM and oligodendrocytes (OLI) partitions were largely similar across tissue types.However, the compositions of MIC, oligodendrocyte precursor cell (OPC), and astrocyte (AST) partitions were unique to certain tissue types.Specifically, the OPC and AST compositions were more similar in pCTX and LGN compared to pWM.On the other hand, the MIC composition was more similar in pWM and pCTX compared to LGN.

Fig
Fig S2 -related to Fig1

Fig
Fig S3 -related to Fig1 also made available for use under a CC0 license.

FigC
Fig S5 -related to Fig2 (B) Stacked bar plots show the relative proportions of transcriptomic ME at the spot resolution across different pathological states (top), as well as the levels of myelin and inflammation (bottom).(C) Dot plots depict the change in spot proportion across subregional white matter (WM) of Samples 5-16.

(
C) Violin plot showing the expression of *CFB (human homolog of marmoset ENSCJAG00000048204), CDKN2A, CYR61, IL16, and VCAM1 across vascular cells in control and EAE.(D) SB/NFR-stained tissue across 16 ROI labeled by the distribution of ME18 (yellow dots) and annotated by lesion age, dated by longitudinal MRI (Methods).
105 and is also made available for use under a CC0 license.(whichwas not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25 Fig S8 -related to Fig2L1 and L2 gene list used to infer cell types

Fig
Fig S9 -related to Fig3 105 and is also made available for use under a CC0 license.(which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25.559371doi: bioRxiv preprint

Fig
Fig S10 -related to Fig3 105 and is also made available for use under a CC0 license.(which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25.559371doi: bioRxiv preprint cell classes and subclusters across tissue type and condition Relative proportion of subclusters across tissue type in EAE Profile of transcription factors enriched in EAE AST MIC OLI OPC GO:0001221_transcription coregulator binding GO:0031962_nuclear mineralocorticoid receptor binding GO:0001046_core promoter sequence−specific DNA binding GO:0001217_DNA−binding transcription repressor activity GO:0031347_regulation of defense response GO:0070721_ISGF3 complex GO:0090077_foam cell differentiation GO:0010742_macrophage derived foam cell differentiation GO:0045637_regulation of myeloid cell differentiation KEGG (C) Venn diagrams illustrate the similarity and diversity of transcription factors that are significantly enriched in EAE compared to control animals across different tissue types for each glial cell class.The elevated transcription factors shared by all tissue types in EAE animals are listed for each cell class, and the shared transcription factors across different cell classes are color-coded accordingly.(D) Dot plot shows selective GO terms enriched for each list of shared transcription factors per glial cell class listed in (C). 105

FigS12
FigS12.related to Fig3.Communication preferences of cells within the white matter (WM) inferred by CellChat.

Fig
Fig S13 -related to Fig4

FigS13
FigS13.related to Fig4.Chord diagrams and network plots summarize the intercellular communication profile of selected pathways. 105
105 and is also made available for use under a CC0 license.(whichwas not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25.559371doi: bioRxiv preprint FigS5.related to Fig2.Matched brain regions of interest across pathological states are color-coded based on different microenvironment (ME) phenotypes.
• • i.GM ii.WM iii.T2 iv.BB v.RP105 and is also made available for use under a CC0 license.(whichwas not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25.559371doi: bioRxiv preprint FigS6.related to Fig2.Spatial organization of transcriptomes within different microenvironments and functional enrichment of gene modules within subregions.
GO:MF HP KEGG 105 and is also made available for use under a CC0 license.(whichwas not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this 105 and is also made available for use under a CC0 license.(whichwas not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC

D Fig S7.continued -related to Fig2
105 and is also made available for use under a CC0 license.(whichwas not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC FigS11.related to Fig3.Comparative analysis of cell types and changes in transcription factors across different tissue types in response to EAE. ( and is also made available for use under a CC0 license.(whichwas not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC ( 105 and is also made available for use under a CC0 license.(whichwas not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC ( which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC also made available for use under a CC0 license.(whichwas not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25.559371doi: bioRxiv preprint ( and is also made available for use under a CC0 license.(which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25.559371doi: bioRxiv preprint 105 and is also made available for use under a CC0 license.(which was not certified by peer review) is the author/funder.This article is a US Government work.It is not subject to copyright under 17 USC The copyright holder for this preprint this version posted September 27, 2023.; https://doi.org/10.1101/2023.09.25.559371doi: bioRxiv preprint