{"paper_id":"28e157bf-d84d-4d97-a919-726ba219c34e","body_text":"Strain-Specific Tropism and Transcriptional Responses of Enterovirus D68 Infection in \nHuman Spinal Cord Organoids \n \nNathânia Dábilla 1¶, Sarah Maya 2¶, Colton McNinch 3, Taylor Eddens 2, Patrick T. Dolan 1*, \nMegan Culler Freeman2* \n \n1Quantitative Virology and Evolution Unit, Laboratory of Viral Diseases, NIH-NIAID Division \nof Intramural Research, Bethesda, MD, USA \n2University of Pittsburgh School of Medicine, Department of Pediatrics, Division of \nInfectious Diseases, Pittsburgh, PA, USA \n \n3Bioinformatics and Computational Bioscience Branch, National Institute of Allergy and \nInfectious Diseases, National Institutes of Health, Rockville, MD, USA \n \n \n \n \n \n \n \n \n \n \n \n$Co-corresponding Authors \nEmails: \nmegan.freeman@pitt.edu (MCF) \npatrick.dolan@nih.gov (PTD) \n¶*These authors contributed equally to this work.  \n \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nAbstract: \nThe mechanisms by which Enterovirus D-68 (EV-D68) infection leads to acute flaccid \nmyelitis (AFM), a severe neurological condition characterized by sudden muscle weakness \nand paralysis, remain poorly understood. To investigate the cellular tropism and infection \ndynamics of EV-D68, we pro filed naive and EV-D68-infected human spinal cord organoids \n(hSCOs) derived from induced pluripotent stem cells (iPSCs) using single-cell RNA \nsequencing (scRNA-seq).  Examining the cellular composition of healthy hSCOs, we found \nthat hSCOs comprise diverse cell types, including neurons, astrocytes, oligodendrocyte \nprogenitor cells (OPCs), and multipotent glial progenitor cells (mGPCs). Upon infection with \ntwo EV-D68 strains, US/IL/14-18952 (a B2 strain) and US/MA/18-23089 (a B3 strain), we \nobserved distinct viral tropism and host transcriptional responses. Notably, US/IL/14-18952 \nshowed a signi\nficant preference for neurons, while US/MA/18-23089 exhibited higher rates \nof infection in cycling astrocytes and OPCs. These findings provide novel insights into the \nhost cell tropism of EV-D68 in the spinal cord, offering insight into the potential mechanisms \nunderlying AFM pathogenesis. Understanding the dynamics of infection at single-cell \nresolution will inform future therapeutic strategies aimed at mitigating the neurological \nimpact of enteroviral infections.  \n2 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nIntroduction \nAcute flaccid myelitis (AFM) is a polio-like illness characterized by muscle weakness and \nparalysis, primarily affecting children [1,2]. Increased cases of AFM suspected to be due to \nenteroviral infection were first recorded in 2014 [3–5] and many of such cases have been \nassociated with Enterovirus D68 (Taxonomy: Enterovirus deconjuncti), or EV-D68 [6–9]. \nSubsequently, EV-D68 and AFM had coinciding biennial outbreaks from 2014 to 2018 [10], \nwith the 2018 AFM outbreak associated with nearly twice as many con firmed AFM cases \ncompared to 2014 [2]. Another spike of cases was expected in 2020, but transmission was \nlikely impeded by isolation policies during the SARS-CoV-2 pandemic [11]. EV-D68 had an \nadditional outbreak in 2022, but AFM cases did not increase as expected [2,12].  \nThe mechanism by which enteroviral infection contributes to the development of AFM is \nunknown. Previous studies have suggested both direct damage to spinal cord neurons after \nviral infection, and subsequent cytotoxic T-cell responses to infected neurons both \ncontribute to disease [9,13]. One fundamental question relevant to EV-D68 pathogenesis is \nthe cell types that contribute to virus replication and production in the CNS. Studies in \nmultiple model systems have demonstrated that EV-D68 can target and replicate in neurons \n[14–18]\n. Astrocyte infection has also been identi fied during EV-D68 infection of murine brain \nslice cultures and primary human astrocytes in vitro [19,20].  \nWe have previously shown that contemporary strains of EV-D68 can replicate in an induced \npluripotent stem cell (iPSC)-derived human spinal cord organoid (hSCO) model, which \nprovides a valuable human-derived, multicellular model in which to explore EV-D68 \npathogenesis [21]. Analysis of marker gene expression suggests hSCOs comprise multiple \ncell types, including neurons and glial cells, but the speci fic cell types present, and which \ncontribute to enteroviral infection in the hSCO model are unknown [21].  \nTo better de fine the cell types infected by EV-D68, and to identify potential strain-speci fic \ndifferences in cell tropism and pathogenesis, we infected hSCOs with two contemporary \nstrains of EV-D68 associated with AFM, US/IL/14-18952 (18952) and US/MA/18-23089 \n(23089). These strains are genetically distinct (18952 is a B2 strain and 23089 is a B3 strain) \nand represent the circulating strains from the 2014 and 2018 outbreak years [22–24] (Fig 1). \nUsing single-cell RNA sequencing (scRNAseq), we captured host and viral transcripts within \ndi\nfferent cell types and subtypes. Our analysis revealed the complex cellular composition of \nhSCOs, which includes neuronal and glial cell lineages. Analyzing viral transcript abundance \nwithin these cell populations demonstrated that the EV-D68 strains exhibit markedly \ndi\nfferent tropisms. Although EV-D68 18952 was primarily associated with neuronal infection, \n23089 exhibited a preference for cycling astrocytes. Subsequent analysis of host cell \ntranscriptional responses in both infected and bystander cell populations revealed further \ndi\nfferences between these strains. Together, these findings clarify the shifting cell tropism of \nEV-D68 strains and provide insight into the mechanism of AFM pathogenesis in a highly \nrelevant human model system.  \nResults \n3 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \n24-days-old spinal cord organoids comprise diverse cell lineages \nPrevious characterization of marker gene and protein expression in hSCOs suggested the \npresence of several fully-di fferentiated cell types, including neuronal lineages, such as \nmotor neurons and interneurons, and roof plate-like structures and neuroepithelium [21,25]. \nHowever, such analyses of gene expression may not identify minor cell types or \ndevelopmental intermediates of speci fic cell lineages, and do not provide quantitative \nmeasures of relative cell abundance. Therefore, to better characterize the diversity and \nabundance of individual cell types in these organoids, we performed scRNAseq. We \ngenerated single-cell suspensions by dissociating pools of 12 hSCOs di\nfferentiated for 24 \ndays, capturing between 8,000 and 10,000 cells per sample for sequencing (Fig 2A).   \nFollowing scRNAseq, we used gene expression pro files from individual cells to perform cell \ntype identi fication. We based our assignments on curated cell types from a previously \npublished study of the developing human spinal cord [26] through “label transfer” [27] (Fig \n2B-C). This analysis revealed that hSCOs exhibit a diverse cellular composition, \nencompassing neuronal lineages, along with astrocytes, oligodendrocytes and their \nprogenitors (OPCs), glial populations such as midplate cells and multipotent glial progenitor \ncells (mGPCs), and vascular leptomeningeal cells (VLMCs), which contribute to blood-brain \nbarrier integrity. \nAstrocytes and cycling astrocytes were identi fied as the majority cell type in the 24-old day \nhSCOs (54.3%). We con firmed their cell identity by assessing expression of known \nastrocyte-specific markers, including SOX9, FGFR3 [28,29] and TOP2A [30] (S1 Fig). \nNeurons comprised the next most common fully-di fferentiated cell type (13.9%), which we \nconfirmed by examining expression of MAP2 and TUBB3 [31,32] (S1 Fig). In addition to \ndifferentiated cell types, mGPCs also made up a considerable proportion of the hSCO \ncomposition (20.1%) (Fig 2B-C, S2-A Fig). These proportions were consistent across \nreplicate pools of hSCOs (S2-B Fig) of 24-days-old hSCOs. Fewer than 10% of cells were \nnot classi\nfied as a speci fic cell type due to low predicted cell type score (Fig 2, \n“Unclassified”). These were not considered for subsequent analysis.  \nFocused reanalysis of the mGPC subset revealed five distinct clusters which may represent \nintermediate states in the di fferentiating organoid (S3-A Fig). Consistent with this \ninterpretation, we have found the proportion of mGPCs are reduced in hSCOs in later \ndevelopment days (data not shown). mGPCs in clusters 0 and 4 display features of neural \nprogenitors and early neuronal di\nfferentiation, with Cluster 0 exhibiting expression of genes \nlinked to mature neuronal identity, while Cluster 4 shows expression patterns indicative of \nactive neurogenesis. Among the expression di fferences, we identified DCX and NEUROG1 \nin Cluster 0, and NKX1-1 and GATA2 in Cluster 4, re flecting their roles in early neural \ndevelopment [33–36]. Cells in Cluster 0 showed high expression of NEUROD4 and ELAVL3 \n[37–39], markers of mature neurons. In contrast, Cluster 4 displays a broader range of \nfunctions, including neurotransmitter synthesis (GAD2 and SLC32A1) [40,41], and cell \nadhesion and signaling (GPR83 and CNTNAP5)  [42–44], indicating a more diverse cell \npopulation (S3-B Fig). Clusters 1, 2, and 3 are likely astrocyte progenitor cells, with diverse \n4 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \ngene expression pro files related to stress response, lipid metabolism, cell signaling, \nextracellular matrix remodeling, DNA repair, cell cycle regulation, neural development, and \ncell polarity. Cluster 2 exhibited markers of proliferation, more consistent with cycling \nastrocytes. \nTwo contemporary EV-D68 strains show distinct tropism in hSCOs \nUpon infection with two distinct strains of EV-D68 (US/IL/14-18952 and US/MA/18-23089), \nwe observed similar patterns of broad cellular susceptibility across the three main cell types \nidentified within the hSCOs (Fig 3-A). Overall, we identi fied 144 viral RNA (vRNA)-positive \ncells out of 8,698 total cells in the 18952-infected organoids, representing 1.65% of the \ncaptured cell population. For strain 23089, we observed 205 vRNA-positive cells out of \n15,805 total cells, corresponding to 1.3% of the captured cells. Although the proportion of \ninfected cells is similar, the two strains exhibited distinct preferences for speci fic cell types. \nEV-D68-18952-infected cells were signi ficantly enriched in neurons (based on permutation \ntests), whereas EV-D68-23089-infected cells were signi ficantly enriched among cycling \nastrocytes and oligodendrocyte progenitor cells (OPC Oligo.) (Fig 3-B). Notably, these \nenrichments were also re flected in the number of viral RNA reads (vRNA) originating from \nthese cell types. Although, of note, we observed significant heterogeneity in vRNA reads per \nindividual cell (Fig 3-C). Together these observations suggest marked strain di fferences in \nhost cell preference in the spinal cord, which may contribute to di fferential pathogenesis \nalthough the consequences are not clear from this observation alone. \nStrain-Specific Enrichment Analysis Reveals Distinct Responses to EV-D68 Infection  \nTo address the distinct cellular and molecular responses triggered by two strains of EV-D68 \n(18952 and 23089), we performed a comprehensive enrichment pathway analysis on \ninfected spinal cord organoids to identify strain-speci fic transcriptional programs. Due to \nthe low proportion of vRNA-positive cells, we combined all cell types to compare gene \nexpression differences based on the infecting strain. While this approach does not resolve \ncell type-speci fic changes, it allows for a more statistically robust comparison of \nstrain-specific di fferences in host response. The normalized enrichment scores (NES) of \nsignificantly enriched pathways (adjusted p < 0.05) revealed three major functional clusters \n(Fig. 4A), capturing pathways where genes are primarily up-regulated (Cluster 1) or \ndown-regulated (Cluster 3) in response to infection with either strain, or pathways primarily \nupregulated in response to 23089 infection relative to 18952-infected (Cluster 2). Notably, \nour analysis did not identify any pathways enriched only in the context of 18952 infection. \nThe most striking di fference between the two strains is in Cluster 2, which corresponded to  \nproliferative and biosynthetic transcriptional programs uniquely enriched in 23089-infected \ncells. This included pathways such as DNA replication, chromosome organization, cell \ncycle, and ribosome biogenesis (Fig. 4A). High NES values were observed for \nDNA-templated DNA replication (NES = 2.07), DNA repair (NES = 1.56), and cell cycle \nprocess (NES = 1.52). Despite 18952-infected cells showing slightly higher expression of \ncell cycle markers per cell (Fig. 4B), the enrichment in 23089 is likely driven by the larger \n5 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nnumber of infected cycling astrocytes, increasing the overall representation of these \nproliferative pathways. 23089-infected cells showed modest upregulation of RNA \nprocessing and splicing pathways across all infected cell types (e.g. SNRNP70, HNRNPA1 \nand SF3B1) (Fig. 4B). These \nfindings suggest that 23089 induces a transcriptionally-active \nand biosynthetically engaged state, possibly re flecting its preferential infection of cycling \nastrocytes. \nClusters 1 and 3 captured similarities in response to the two strains. Pathways identi fied in \nCluster 1 included neuronal development and communication—such as synapse assembly, \naxon development, neurotransmitter secretion, G protein-coupled receptor signaling, and \nneurogenesis—which were positively enriched in both 18952- and 23089-infected cells, \nconsistent with the neurotropic nature of EV-D68. However, strain 18952 showed higher \nenrichment for synaptic signaling (NES = 1.87 vs. 1.13) and neurotransmitter transport (NES \n= 1.95 vs. 1.13), suggesting a more robust neuronal response (Fig 4A). Cluster 3 \nencompassed pathways related to translation, metabolism, and developmental signaling, \nwith 18952 showing stronger negative enrichment, including cytoplasmic translation (NES = \n–2.22), AMP metabolic process (NES = –2.26), and ubiquitin ligase regulation (NES = –2.30). \nWhile genes encoding ribosomal proteins (e.g. RPL10, RPLP1 and RPS6) appeared \ndownregulated in neurons from both infected and mock conditions (Fig 4B), broader \nsuppression of metabolic homeostasis pathways was more pronounced in the infected \ncells. Together, these clusters suggest that both strains impact neuronal function and basal \ncellular processes, but with greater intensity in 18952-infected cells. \nThese transcriptional di fferences align with the distinct cellular tropism observed between \nthe two EV-D68 strains. In organoids infected with 18952, neurons represented the major \nsignificantly enriched infected cell population, which is consistent with the prominent \nupregulation of synaptic signaling, axonal development, and neurotransmitter-related \npathways. The robust enrichment of these neuronal processes suggests a direct viral \nimpact on neurons, potentially contributing to neuronal dysfunction or degeneration. In \ncontrast, 23089 predominantly infected cycling astrocytes and, to a lesser extent, \noligodendrocyte progenitor cells (OPCs). This cell-type speci\nficity is mirrored by the \nenrichment of DNA replication, cell cycle progression, and RNA processing pathways in \n23089-infected cells—hallmarks of transcriptionally active, proliferative glial populations. \nThus, the pathway signatures not only re flect divergent viral-host interactions but also \nunderscore how strain-speci fic cellular targeting shapes the overall transcriptional \nlandscape of infected organoids.  \nDistinct Transcriptional Impacts of EV-D68 Strains Across Bystander Cell Types \nDifferential gene expression (DEG) analysis was performed on bystander cells strati fied by \ncell type. A summary of the number of signi ficant DEGs identified per cell type is presented \nin (Fig. 5A). This panel re flects DEGs filtered by both adjusted p-value < 0.05 and absolute \nlog2 fold-change greater than 0.25, ensuring that the displayed genes represent biologically \nrelevant transcriptional shifts. Signi ficant DEG counts were observed for neurons, \n6 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nastrocytes, and cycling astrocytes, whereas OPCs did not show any DEGs that met these \ncriteria. Accordingly, they were analyzed similarly to the vRNA-positive cells.  \nThe transcriptional response in bystanders varied notably according to both cell type and \nviral strain. Among the di fferent cell types, neurons demonstrated a short list of DEGs, while \nastrocytes and especially cycling astrocytes showed a more pronounced transcriptional \nresponse. Neurons and astrocytes displayed a greater number of DEGs in response to the \n23089 strain, suggesting a more pronounced impact of this viral strain on these cell \npopulations. In contrast, cycling astrocytes exhibited more DEGs associated with the 18952 \nstrain. However, the magnitude of transcriptional changes di\nffered between strains: DEGs \nassociated with 18952 tended to have more modest fold changes, whereas 23089 DEGs \nspanned a broader range of fold changes (Fig 5D). These patterns indicate that both the \nextent and intensity of the bystander transcriptional response di ffer depending on cell type \nand viral strain. \nNeurons displayed a comparatively limited bystander transcriptional response, with fewer \ndifferentially expressed genes than astrocytes or cycling astrocytes. In the 18952 condition, \nupregulated genes included SLC5A7, essential for acetylcholine synthesis and synaptic \ntransmission [45], which is consistent with the cholinergic nature of spinal motor neurons. \nECEL1, implicated in motor neuron axon development [46], was also elevated. Other \ntranscripts, such as HOXA4 and CCBE1, related to positional identity and extracellular \nmatrix organization [47,48], suggest a role in preserving neuronal identity and structure. \nIn contrast, 23089-exposed neurons upregulated stress-responsive and protein quality \ncontrol genes. These included PTPN11 [49], CHORDC1, and PRKCSH, associated with \nMAPK signaling, protein folding, and ER stress [50,51]. Increases in MT-ND6 and PSMC4 \nexpression, involved in mitochondrial respiration and proteasomal degradation [52,53], may \nreflect compensatory responses to maintain homeostasis (Fig. 5B). These trends suggest \nthat while 18952 maintains neuronal identity, 23089 induces a mild stress-adaptive \ntranscriptional state.  \nAstrocytes exhibited a broader transcriptional response than neurons, with distinct gene \nexpression patterns between strains. In the 18952 condition, upregulated genes included \nCHST9 and PDGFA, involved in extracellular matrix remodeling and glial signaling, \nalongside HRH3 and SLC2A3, linked to histaminergic modulation and glucose transport. \nAdditional genes such as BCL3, PFKFB4, and MT1X suggested mild activation of redox \nregulation and metabolic homeostasis. Together, these changes re flect a modestly reactive \nbut structurally supportive astrocyte state. In contrast, astrocytes exposed to the \n23089-strain upregulated a distinct set of genes, including HSPA1A, CNTNAP2, PCDH19, \nand ROBO2, which are involved in synaptic regulation, axon-glia interaction, and cellular \nstress responses. The transcriptional pro file also included ERBB4, LINGO2, and TOP2A, \nfurther implicating altered signaling and adhesion programs. These di fferences suggest that \n23089-exposed astrocytes enter a more transcriptionally active state, potentially a ffecting \n7 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \ncommunication with the neuronal environment. The overall transcriptional divergence \nreinforces the strain-specific modulation of astrocyte identity in the bystander context (Fig. \n5C). \nCycling astrocytes displayed a unique transcriptional pro file distinct from both neurons and \nnon-cycling astrocytes, marked by broader identity and functional shifts between strains. In \nthe 18952 condition, upregulated genes included SERPINB9, DPAGT1, and LOX, along with \nhistone-associated transcripts such as H2AC13, H4C3, and HIST1H1B. These changes \nsupport a transcriptional program favoring immune regulation, extracellular matrix \norganization, and cell cycling. 23089-exposed cycling astrocytes showed a marked shift in \nidentity, with upregulation of neuronal and stress-related genes including STMN2, NEFM, \nDCX, and CRABP1. Additional activation of HSPA1A, JUN, and ERBB4 suggested elevated \nstress signaling and cytoskeletal remodeling, consistent with transcriptional dysregulation. \nThe emergence of neuronal lineage markers in this glial progenitor population may re\nflect \nstress-induced unexpected activation of developmental genes not typically active in \nastrocyte progenitors [54].  \nNotably, this mixed transcriptional phenotype in cells predicted to be cycling astrocytes  \nmay be linked to the fact that these cells are also the primary targets of infection by the \n23089 strain. Even in bystander cells, signaling cues from nearby infected cells could \nperturb their transcriptional program and disrupt normal progenitor dynamics. Compared to \nnon-cycling astrocytes, which maintained more canonical glial features, the 23089-exposed \ncycling astrocytes exhibited broader deviations in transcription. These distinctions highlight \nthe unique vulnerability of proliferative astrocyte populations to strain-speci\nfic viral influence \n(Fig. 5D). \nIn OPC Oligodendrocytes, pathway enrichment analysis revealed distinct transcriptional \nprograms between strains, even in the bystander population. Due to the relatively low \nnumber of OPC Oligo. in this group, we adopted the same strategy used for infected cell \nanalysis, performing enrichment pathways analysis to enable broader biological \ninterpretation since no gene by itself stands out from our DEG analysis. OPCs exposed to \nthe 18952-strain showed upregulation of pathways associated with lipid metabolism and \nextracellular matrix organization, including foam cell di\nfferentiation, regulation of cholesterol \nefflux, and aminoglycan biosynthetic process. These processes are important for membrane \ndynamics and maintenance of oligodendrocyte identity, suggesting that in the 18952 \ncondition they may retain a more metabolically stable and structurally supportive state. \nIn contrast, OPCs exposed to the 23089-strain were enriched for pathways related to cell \ncycle regulation (mitotic spindle elongation, cytokinesis), microtubule organization, and \ndevelopmental signaling, including glial cell fate speci fication. These transcriptional \nsignatures are consistent with a disturbed cellular state, which could potentially re flect \ndirect infection-related stress or dysregulation, as we observe for cycling astrocytes. This is \n8 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nin line with the fact that OPCs, alongside cycling astrocytes, represent a signi ficant infected \npopulation in the 23089 condition. The enrichment of mitotic and developmental signaling \nprograms may indicate that viral presence disrupts the normal proliferative or \nlineage-committed state of these glial progenitors. Together, these findings suggest that \nOPCs respond to 23089 exposure with transcriptional shifts that may compromise their \nhomeostatic roles (Fig. 5E). \nFinally, we noted a lack of altered expression of innate immune genes, particularly interferon \nstimulated genes (ISGs), in our organoid model. Further analysis revealed we could detect \nexpression of genes known to be regulated by Type I, II, and III interferon in both infected \nand uninfected organoids suggesting the lack of induction was not due to limited detection \nin our sequencing (Supp. Fig. 4). Although, more robust immune responses may occur at \nlater time points in infection (beyond 2dpi), this lack of ISG induction is consistent with \nprevious studies highlighting the e\nffectiveness of viral innate immune suppression in \nEV-D68 infection, through protease cleavage of innate immune signalling factors and the \naction of other viral proteins [55,56].  \nDiscussion \nOur analysis of 24-day old organoids showed a diverse cell composition, with both neuronal \nand glial lineages, consistent with previous IF marker analysis [25]. While we identi fied \nneurons, mGPCs, OPCs, VLMCs, and more, the most common cell type in the hSCO at this \nage was cycling astrocytes and astrocytes (54.3%). Given the relatively nascent nature of \nthe hSCO model, this majority may be due to astrocytes' vital role in developing and \nmaintaining neuronal functions throughout spinal cord development [57–60]. mGPCs \n(20.1%) also have the capacity to differentiate into both astrocytes and oligodendrocytes so \ncellular composition in more mature hSCO may shift [61]. Although this cell population has \nnot yet committed to a di fferential pathway, they are still more specialized and committed \nthan a broad progenitor cell. 24-day old organoids allow for the growth and development of \nCNS cell types without compromising on overall cell viability [21]. The abundance of CNS \ncell types present and interacting in hSCO allow for us to better understand EV-D68 tropism \nand dynamics in the complex human spinal cord.  \n \nTo understand which CNS cell-types are infected by EV-D68, we performed scRNA-seq on \nhSCO infected with contemporary EV-D68 strains US/IL/14-18952 and US/MA/18-23089. \nWe found that while both strains had the capacity to infect neuronal and glial cell lineages, \n18952 preferentially infected neurons while 23809 preferentially infected cycling astrocytes \nand to a lesser extent, oligodendrocyte precursor cells. These di\nfferences extended to \nbystander cells — i.e., uninfected cells within infected hSCOs — where we observed more \nDEGs in 18952-bystander cycling astrocytes, and even larger fold changes in \n23089-bystander neurons and astrocytes. These \nfindings indicate that glial cell populations \nsuch as astrocytes and oligodendrocytes play an important role in EV-D68 associated AFM \npathogenesis, varying amongst di fferent strains and clade classi fication. Previous studies \n9 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nhave implicated spinal cord neurons in EV-D68-induced paralysis, but our results suggest \nadditional cell types are infected and transcriptionally altered, potentially playing \nunidentified roles in EV-D68 pathogenesis [13–17].  \n \nAstrocytes participate in both innate and adaptive immune system, such as regulating the \nrelease of cytokines and chemokines, thereby leading to antigen presentation that lead to \nrecruitment of helper T-cells in pathogen or damage-a ffected areas of the CNS [62–66]. \nSeveral of the di fferentially expressed genes identi fied in bystander astrocytes and cycling \nastrocytes, particularly in the 23089 condition, are associated with stress signaling \npathways. These changes, combined with the disruption of glial identity in 23089-exposed \ncycling astrocytes, may create an environment more permissive to immune cell recruitment. \nIn contrast, the transcriptional pro file in 18952-exposed astrocytes suggested a more \nlimited activation of broad immune-modulatory programs, although cycling astrocytes did \nexpress immunoregulatory markers such as CD70 and SERPINB9, potentially re\nflecting a \nmore localized or cell-type stress response rather than a widespread in flammatory \nactivation. These findings raise the possibility that di fferential modulation of innate immune \nsignaling by each strain could further shape their pathogenic outcomes in vivo. Since we do \nnot know yet to what degree the immune response mediates EV-D68-associated AFM, \nfurther studies are needed to assess astrocyte function and subsequent immune response \nduring EV-D68 infection. Additionally, because oligodendrocytes form myelin to aid neuron \nconductivity and communication, its functionality upon EV-D68 infection should also be \nassessed [67,68]. \n \nWe also identified strain-specific differences in cellular tropism within hSCO, indicating that \nthere has been a change in EV-D68 viral infection dynamics from 2014 to 2018. \nEV-D68-18952 strain primarily infected neurons, while 23089 preferentially infected cycling \nastrocytes and, to a lesser extent, OPCs. Together, these results suggest that the two \nstrains may induce distinct forms of cellular vulnerability — 18952 through direct neuronal \ntargeting and 23089 through glial destabilization. While further studies will be needed to \ndetermine the long-term impact of these responses, the transcriptional divergence \nobserved across cell types points to fundamentally di\nfferent modes of pathogenesis. These  \ndifferences may be broadly attributed to clade speci fic pathogenic di fferences between \n18952 (B2) and 23089 (B3) or be pathogenic variations between two speci fic isolates. While \nboth viruses are expected to cause the same clinical paralysis phenotype, the differences in \ntropism may have caused varying mechanisms of pathogenesis and we are not able to \nassess clinical di\nfferences between these two isolates.  \n \nAdditional strain specific differences are reflected in the transcriptional impacts in bystander \ncells. Astrocytes and cycling astrocyte bystander cells had significantly more transcriptional \nchanges compared to neuronal bystanders, with 23089 having more impact on astrocyte \nand neuron bystanders and 18952 having more of an impact on cycling astrocyte \nbystanders. Many of the upregulated genes in 23089 bystander cells were cellular stress \nindicators, primarily found in astrocytes. These cellular responses may also indicate a \n10 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \ndownstream loss in communication or homeostasis in neurons, as astrocytes are integral to \nmaintaining neuronal function and stability. 23089 had a larger fold transcriptional shift in \nbystander cells compared to 18952. Although these are only two isolates, these di fferences \nmay point to clade- or outbreak-speci fic variation in EV-D68 pathogenesis. The higher \nnumber of AFM cases in 2018 compared to 2014 could re flect such viral di fferences, but \nmay also stem from multifactorial causes, including co-circulation of other enteroviruses \n(e.g., EV-A71) or improved case detection and reporting [69].  \n \nThe OPC bystander population did not have any signi ficant shifts, 18952 exposed OPCs \nwere more stable and 23089 exposed OPCs had evidence of a stress response that \ninterfered with glial cell di fferentiation. EV-D68 may be interfering with oligodendrocyte \ndifferentiation and maturation in the hSCO and thus prevent a fully functioning \noligodendrocyte population. Considering the hSCO population of oligodendrocytes were \nprecursor cells and still differentiating at the time of EV-D68 infection, cellular response may \ndiffer in mature oligodendrocytes. To better understand oligodendrocyte’s response to \nEV-D68 infection, studies in more aged hSCO will be necessary.   \n \nWhile hSCO scRNAseq has allowed us to further understand EV-D68 infection dynamics in \nthe spinal cord, it does have its limitations. Since hSCO cells are relatively immature and \nrepresent a developing human spinal cord, EV-D68 infection dynamics may di\nffer in hSCO \nto that of a child. Despite this limitation, its multicellular complexity and physiological \nrelevance for AFM has and can let us learn more about EV-D68 pathogenesis in the CNS. \nThis limitation in hSCO also lends itself to be an advantage when looking at enteroviruses \nthat target neonates in order to better understand their tropism and pathogenesis.  \n \nAnother limitation in our study is the relatively low number of infected cells identi fied for \nboth EV-D68 strains. Despite leveraging an aggregated, bulk analysis approach to enhance \nstatistical power, the small fraction of virus-positive cells could potentially limit the detection \nof more subtle transcriptional changes and low-abundance cell populations responding to \ninfection. This may be particularly relevant for less abundant cell types, such as OPCs. \nMoreover, oligodendrocytes are known to be particularly fragile during tissue dissociation, \npotentially leading to underrepresentation in the \nfinal dataset. Additionally, the low number \nof infected cells observed could also re flect loss during dissociation or capture, particularly \nif infected cells were damaged and lysed during infection or processing.  \n \nDespite these limitations, our study provides valuable insights into the cellular and \nmolecular landscape of EV-D68 infection in human spinal cord organoids, revealing distinct \ncell-type tropism for two contemporary strains. By integrating pathway enrichment analysis \nwith single-cell transcriptomic pro filing, we demonstrate that EV-D68-18952 exhibits a \npronounced preference for neurons, driving synaptic signaling and axonal development \npathways, while EV-D68-23089 predominantly targets cycling astrocytes, triggering \ntranscriptional programs associated with cell cycle progression and RNA processing. These \nfindings represent a step forward in understanding the strain-specific interactions of EV-D68 \n11 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nwith neural populations, which may have implications for viral spread and \nneuropathogenesis. Furthermore, our use of spinal cord organoids as a model provides a \nphysiologically relevant system to dissect host–virus interactions at single-cell resolution, \nunderscoring the utility of this platform for studying neurotropic viruses.  \n12 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nMaterials and Methods \nViruses and cells \nEV-D68 US/IL/14-18952 (CDC) and US/MA/18-23089 (CDC) strains were propagated using \nHeLa cells incubated at 33°C and 5% CO 2 and purified using sucrose-cushion as previously \ndescribed [70]. These stocks were previously sequenced to con firm their identity with VP1 \nprimers3.  \nHeLa 7b (ATCC, CCl-2) cells were maintained in MEM medium (ThermoFisher, 11095-072 ), \nsupplemented to contain 5% FBS (Phenonix Scienti fic, PS-100), 1% penicillin/streptomycin \n(Corning, 30-002-Cl), and 1% NeAA (Corning, 25-025-Cl). Cells were grown at 37℃ and 5% \nCO2. \nHuman iPSC line SCTi003A  (STEMCELL Technologies, 200-0511) was maintained in \nmTeSRTM Plus medium (STEMCELL Technologies, 100-0276), supplemented with 10 µM \nY-27632 (Tocris, 1254). They were seeded and passaged in flasks coated with 150 µg/mL \nCultrex (R&D Systems, 3434-005-02). The 3-DiSC hSCO were propagated and \ndifferentiated as described in Aguglia et al. [21] for up to 24 days.  \nhSCO infections and dissociations \nhSCOs were infected in pools of 12 organoids each and inoculated with virus at 10 5 \nPFU/pool. After 1 hour of incubation at room temperature, hSCOs were washed 3X with \nPBS and moved to new wells before incubation with fresh medium. No further media \nchanges were performed for the rest of the experiment. The EV-D68 infected pool was \nincubated at 33 oC for 48 hours post infection (hpi). After 48hpi, the pools were dissociated \nto a single-cell suspension with Accumax (Invitrogen). They were incubated for 15 minutes \nin the water bath at 37℃ , gently mixed, and then proceeded with the proposed 10X \nGenomics protocol (Cell Preparation Guide - CG00053 Rev C). Once the cell's \nconcentration was achieved the cells were moved to ice. \nSingle cell RNAseq cDNA library generation \nAll samples were calculated to achieve ~5000-8000 targeting cells in the single-cell \npreparations using Chromium Next GEM Single Cell 5’ standard kit. For the preparation of \nthe cDNA and sequencing library generation, we followed the instructions from the user \nguide Chromium Next GEM Single Cell 5’ Reagent Kit v2 (Dual index). All other steps were \nfollowed to produce cDNA and subsequent Illumina sequencing library for single cell \nsequencing. The illumina library preparation was submitted to quality control in the \nTapeStation D1000 high sensitivity for size distribution and DNA concentration was \nmeasured by Qubit High Sensitivity dsDNA kit. The molar concentration of the libraries were \ndetermined and the samples were diluted for sequencing according to illumina sequencing \nprotocol. We aimed to sequence each library to achieve ~50,000 reads per cell. \nCellRanger \n13 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nCellRanger v7.0.0 was used [71] to align reads against a composite genome, which \nencompassed both the GRCh38 human reference genome and the Enterovirus genome \ncorresponding to the speci fic strain identity of the sample. Feature-barcode matrices were \ngenerated using GENCODE v44 GRCh38 gene models and the viral strain’s genome as a \nsingle ORF. Default parameters of the ‘cellranger count’ function were used. \nSeurat \nSeurat v5.0.2 was used [72]. CellRanger gene counts were made compatible across \ndifferent strains by setting the name of the gene encompassing all viral reads for each \nsample to ‘Viral-Gene’. Cells with mitochondrial gene expression > 25% or detected gene \ncounts < 1000 were removed. Gene counts of filtered cells were normalized using \nSCTransform, while specifying the ‘vars.to.regress’ parameter to the ‘S.Scores’ and \n‘G2M.Scores’ obtained from the ‘CellCycleScoring’ function, utilizing Seurat’s ‘cc.genes’ \ncell cycle gene list. Reciprocal PCA integration analyses were then performed to generate \nintegrated datasets for EV-D68. Dimensionality reduction of the integrated dataset was \nperformed using the ‘RunPCA’ function and ‘RunUMAP’ function using the top 30 principal \ncomponents. Cell cluster analysis was performed using the ‘FindNeighbors’ function using \nthe top 30 principal components and the ‘FindClusters’ function using a ‘resolution’ of 0.1. \nCell type annotations were predicted using cell label transfer with Seurat’s \nFindTransferAnchors, TransferData, and AddMetaData functions using a previously \nannotated spinal cord scRNAseq dataset [26] as our reference, considering “true” cell types \nabove 0.5 prediction.score.max. The cells that were below this threshold were assigned as \nNA or Unclassified.  \nAssignment of infection status \nViral read percentages were calculated as the fraction of total UMIs per cell mapping to the \nviral gene (i.e. ‘Viral-Gene’) and were computed using the ‘PercentageFeatureSet’ Seurat \nfunction. To classify cells as 'Infected' or 'NonInfected', a Poisson test was applied to these \npercentages using the ‘estimateNonExpressingCells’ function from the SoupX R package \n[73], using an FDR of 0.05. This function accounts for ambient RNA contamination in the \nsample, which was estimated using SoupX's 'autoEstCont' function with 't\nfi\ndfMin' set to 1.0 and \n'soupQuantile' set to 0.9. \nPermutation-based enrichment analysis \nTo assess whether speci fic cell types were signi ficantly enriched in infected (vRNA+ cells) \nversus non-infected (Bystander cells) conditions, we performed a permutation test on cell \ncount distributions across infection states. For each viral strain, we computed contingency \ntables comparing observed cell type frequencies across infection status (FDR < 0.05). To \nestablish a null distribution, we randomly permuted cell type labels 10,000 times while \npreserving the infection status labels and recomputed the contingency tables for each \niteration. Median values and 95% con\nfidence intervals (2.5th and 97.5th percentiles) were \nderived from the permuted distributions. Observed counts exceeding the upper con fidence \n14 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \ninterval were considered signi ficantly enriched (p value < 0.05). Enrichment scores were \ncalculated as the ratio of observed to permuted median counts for each cell type and \ncondition. Analyses were conducted in R using the data.table, Seurat, and SeuratObject \npackages. \nDifferential gene expression analysis \nDifferential gene expression (DEG) analysis was conducted using the FindMarkers function \nin Seurat v5. Two complementary approaches were applied. First, a pseudo-bulk-style \nanalysis was performed by grouping cells according to their combined infection status and \nviral strain (e.g., “18952_Infected” vs. “18952_NotInfected”), enabling comparisons across \naggregated conditions. This approach was used because the number of infected cells was \ninsufficient for robust strati fication by cell type. Second, a cell-type-speci fic DEG analysis \nwas performed using only bystander cells, focusing on neurons, astrocytes, cycling \nastrocytes, and oligodendrocyte lineage cells (OPC Oligo.) to compare transcriptional \nresponses between EV-D68 strains 18952 and 23089. All analyses were conducted using \nthe RNA assay, with normalization via NormalizeData. DEG identification used a minimum \nexpression threshold of 10% (min.pct = 0.1), no log fold-change cuto ff (logfc.threshold = 0), \nand a minimum of 50 cells per group. Fold changes were calculated as log2-transformed \nvalues, and signi ficance was assessed using Bonferroni-adjusted p-values. Genes with \nadjusted p value < 0.05 and | log₂FC | > 0.25 were considered signi ficantly differentially \nexpressed. Visualizations, including volcano plots and DEG count barplots by cell type and \ndirection, were generated using ggplot2, ggrepel, and patchwork. \nPathway enrichment analysis \nGene set enrichment analysis (GSEA) was performed using the fgsea package to identify \nbiological pathways enriched in infected conditions or speci fic cell types, even in the \nabsence of signi ficantly di fferentially expressed genes. For certain DEG comparisons, \nparticularly those involving low-abundance populations (e.g., infected cells or bystander \nOPC Oligo.), no genes met the adjusted p-value threshold for signi ficance. In these cases, \nthe full ranked DEG lists based on average log2 fold change were used as input for GSEA, \nenabling the detection of coordinated pathway-level shifts. Ranked gene lists were \ngenerated from pairwise comparisons of infection conditions (e.g., 18952-Infected vs. \nMock) and bystander cell types (e.g., 23089 vs. 18952 within OPC Oligo.), and enrichment \nwas computed using 100 million permutations for robust estimation. Gene sets were \nsourced from the Gene Ontology Biological Process (GO-BP) category via msigdbr. \nNormalized enrichment scores (NES) were computed for each condition, and pathways with \nadjusted p-values < 0.05 were considered signi ficantly enriched. Heatmaps of NES values \nacross conditions were generated using ComplexHeatmap, and selected pathway \ncomparisons were visualized with ggplot2 and ggrepel. \nData Availability \nData from the scRNAseq analysis is deposited in Gene Expression Omnibus (GEO) \n15 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \n(Accession number GSE292051). \nFunding Statement \nThis work was supported by the Intramural Research Program of the National Institute of \nAllergy and Infectious Diseases at the National Institutes of Health - Project Number \n1ZIAAI001360 (PTD). MCF receives support from National Institutes of Health K08AI171177. \n16 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nReferences \n1.  Murphy OC, Messacar K, Benson L, Bove R, Carpenter JL, Crawford T, et al. Acute \nflaccid myelitis: cause, diagnosis, and management. Lancet. 2021;397: 334–346. \n2.  CDC. AFM Cases & Outbreaks. In: Acute Flaccid Myelitis (AFM) [Internet]. 7 Feb 2025 \n[cited 5 Mar 2025]. Available: \nhttps://www.cdc.gov/acute-flaccid-myelitis/cases/index.html \n3.  Messacar K, Asturias EJ, Hixon AM, Van Leer-Buter C, Niesters HGM, Tyler KL, et al. \nEnterovirus D68 and acute flaccid myelitis—evaluating the evidence for causality. \nLancet Infect Dis. 2018;18: e239–e247. \n4.  Messacar K, Abzug MJ, Dominguez SR. 2014 outbreak of enterovirus D68 in North \nAmerica. J Med Virol. 2016;88: 739–745. \n5.  Aliabadi N, Messacar K, Pastula DM, Robinson CC, Leshem E, Sejvar JJ, et al. \nEnterovirus D68 Infection in Children with Acute Flaccid Myelitis, Colorado, USA, 2014. \nEmerg Infect Dis. 2016;22: 1387–1394. \n6.  McKay SL, Lee AD, Lopez AS, Nix WA, Dooling KL, Keaton AA, et al. Increase in acute \nflaccid myelitis - United States, 2018. MMWR Morb Mortal Wkly Rep. 2018;67: \n1273–1275. \n7.  Schubert RD, Hawes IA, Ramachandran PS, Ramesh A, Crawford ED, Pak JE, et al. \nPan-viral serology implicates enteroviruses in acute flaccid myelitis. Nat Med. 2019;25: \n1748–1752. \n8.  Mishra N, Ng TFF, Marine RL, Jain K, Ng J, Thakkar R, et al. Antibodies to \nenteroviruses in cerebrospinal fluid of patients with acute flaccid myelitis. MBio. \n2019;10. doi:10.1128/mBio.01903-19 \n9.  Vogt MR, Wright PF, Hickey WF, De Buysscher T, Boyd KL, Crowe JE Jr. Enterovirus \nD68 in the Anterior Horn Cells of a Child with Acute Flaccid Myelitis. N Engl J Med. \n2022;386: 2059–2060. \n10.  Shah MM, Perez A, Lively JY, Avadhanula V, Boom JA, Chappell J, et al. Enterovirus \nD68-associated acute respiratory illness ─ New Vaccine Surveillance Network, United \nStates, July-November 2018-2020. MMWR Morb Mortal Wkly Rep. 2021;70: \n1623–1628. \n11.  Olsen SJ, Winn AK, Budd AP, Prill MM, Steel J, Midgley CM, et al. Changes in influenza \nand other respiratory virus activity during the COVID-19 pandemic - United States, \n2020-2021. MMWR Morb Mortal Wkly Rep. 2021;70: 1013–1019. \n12.  New Vaccine Surveillance Network Collaborators, Hall AJ. 2022. Increase in acute \nrespiratory illnesses among children and adolescents associated with rhinoviruses and \nenteroviruses. \n13.  Woods Acevedo MA, Lan J, Maya S, Jones JE, Williams JV, Freeman MC, et al. \nImmune cells promote paralytic disease in mice infected with enterovirus D68. bioRxiv. \n2024. doi:10.1101/2024.10.14.618341 \n17 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \n14.  Hixon AM, Yu G, Leser JS, Yagi S, Clarke P, Chiu CY, et al. A mouse model of paralytic \nmyelitis caused by enterovirus D68. PLoS Pathog. 2017;13: e1006199. \n15.  Brown DM, Hixon AM, Oldfield LM, Zhang Y, Novotny M, Wang W, et al. Contemporary \nCirculating Enterovirus D68 Strains Have Acquired the Capacity for Viral Entry and \nReplication in Human Neuronal Cells. MBio. 2018;9. doi:10.1128/mBio.01954-18 \n16.  Rosenfeld AB, Warren AL, Racaniello VR. Neurotropism of Enterovirus D68 isolates is \nindependent of sialic acid and is not a recently acquired phenotype. MBio. 2019;10. \ndoi:10.1128/mBio.02370-19 \n17.  Feng M, Guo S, Fan S, Zeng X, Zhang Y, Liao Y, et al. The preferential infection of \nastrocytes by Enterovirus 71 plays a key role in the viral neurogenic pathogenesis. Front \nCell Infect Microbiol. 2016;6: 192. \n18.  Wang Y-F, Chou C-T, Lei H-Y, Liu C-C, Wang S-M, Yan J-J, et al. A mouse-adapted \nEnterovirus 71 strain causes neurological disease in mice after oral infection. J Virol. \n2004;78: 7916–7924. \n19.  Rosenfeld AB, Warren AL, Racaniello VR. Neurotropism of Enterovirus D68 Isolates Is \nIndependent of Sialic Acid and Is Not a Recently Acquired Phenotype. mBio. American \nSociety for Microbiology; \n20.  Liu X, Li H, Li Z, Gao D, Zhou J, Ni F, et al. MFSD6 is an entry receptor for respiratory \nenterovirus D68. Cell Host Microbe. 2025;33: 267–278.e4. \n21.  Aguglia G, Coyne CB, Dermody TS, Williams JV, Freeman MC. Contemporary \nenterovirus-D68 isolates infect human spinal cord organoids. MBio. 2023;14: \ne0105823. \n22.  Hadfield J, Megill C, Bell SM, Huddleston J, Potter B, Callender C, et al. Nextstrain: \nreal-time tracking of pathogen evolution. Bioinformatics. 2018;34: 4121–4123. \n23.  Sagulenko P, Puller V, Neher RA. TreeTime: Maximum-likelihood phylodynamic \nanalysis. Virus Evol. 2018;4: vex042. \n24.  Brown BA, Nix WA, Sheth M, Frace M, Oberste MS. Seven strains of Enterovirus D68 \ndetected in the United States during the 2014 severe respiratory disease outbreak. \nGenome Announc. 2014;2. doi:10.1128/genomeA.01201-14 \n25.  Ogura T, Sakaguchi H, Miyamoto S, Takahashi J. Three-dimensional induction of \ndorsal, intermediate and ventral spinal cord tissues from human pluripotent stem cells. \nDevelopment. 2018;145. doi:10.1242/dev.162214 \n26.  Andersen J, Thom N, Shadrach JL, Chen X, Onesto MM, Amin ND, et al. Single-cell \ntranscriptomic landscape of the developing human spinal cord. Nat Neurosci. 2023;26: \n902–914. \n27.  Stuart T, Butler A, Hoffman P, Hafemeister C, Papalexi E, Mauck WM 3rd, et al. \nComprehensive integration of single-cell data. Cell. 2019;177: 1888–1902.e21. \n28.  Stolt CC, Lommes P, Sock E, Chaboissier M-C, Schedl A, Wegner M. The Sox9 \ntranscription factor determines glial fate choice in the developing spinal cord. Genes \n18 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nDev. 2003;17: 1677–1689. \n29.  Pringle NP, Yu W-P, Howell M, Colvin JS, Ornitz DM, Richardson WD. Fgfr3 expression \nby astrocytes and their precursors: evidence that astrocytes and oligodendrocytes \noriginate in distinct neuroepithelial domains. Development. 2003;130: 93–102. \n30.  Lee JH, Berger JM. Cell cycle-dependent control and roles of DNA topoisomerase II. \nGenes (Basel). 2019;10: 859. \n31.  Menezes JR, Luskin MB. Expression of neuron-specific tubulin defines a novel \npopulation in the proliferative layers of the developing telencephalon. J Neurosci. \n1994;14: 5399–5416. \n32.  Dehmelt L, Halpain S. The MAP2/Tau family of microtubule-associated proteins. \nGenome Biol. 2005;6: 204. \n33.  Gleeson JG, Lin PT, Flanagan LA, Walsh CA. Doublecortin is a microtubule-associated \nprotein and is expressed widely by migrating neurons. Neuron. 1999;23: 257–271. \n34.  Ma Q, Fode C, Guillemot F, Anderson DJ. Neurogenin1 and neurogenin2 control two \ndistinct waves of neurogenesis in developing dorsal root ganglia. Genes Dev. 1999;13: \n1717–1728. \n35.  Nardelli J, Thiesson D, Fujiwara Y, Tsai FY, Orkin SH. Expression and genetic interaction \nof transcription factors GATA-2 and GATA-3 during development of the mouse central \nnervous system. Dev Biol. 1999;210: 305–321. \n36.  Schubert FR, Fainsod A, Gruenbaum Y, Gruss P. Expression of the novel murine \nhomeobox gene Sax-1 in the developing nervous system. Mech Dev. 1995;51: 99–114. \n37.  Okano HJ, Darnell RB. A hierarchy of Hu RNA binding proteins in developing and adult \nneurons. J Neurosci. 1997;17: 3024–3037. \n38.  Miyata T, Maeda T, Lee JE. NeuroD is required for differentiation of the granule cells in \nthe cerebellum and hippocampus. Genes Dev. 1999;13: 1647–1652. \n39.  Mulligan MR, Bicknell LS. The molecular genetics of nELAVL in brain development and \ndisease. Eur J Hum Genet. 2023;31: 1209–1217. \n40.  Erlander MG, Tobin AJ. The structural and functional heterogeneity of glutamic acid \ndecarboxylase: a review. Neurochem Res. 1991;16: 215–226. \n41.  McIntire SL, Reimer RJ, Schuske K, Edwards RH, Jorgensen EM. Identification and \ncharacterization of the vesicular GABA transporter. Nature. 1997;389: 870–876. \n42.  Müller TD, Müller A, Yi C-X, Habegger KM, Meyer CW, Gaylinn BD, et al. The orphan \nreceptor Gpr83 regulates systemic energy metabolism via ghrelin-dependent and \nghrelin-independent mechanisms. Nat Commun. 2013;4: 1968. \n43.  Chatterjee M, Schild D, Teunissen CE. Contactins in the central nervous system: role in \nhealth and disease. Neural Regen Res. 2019;14: 206–216. \n44.  Gomes I, Bobeck EN, Margolis EB, Gupta A, Sierra S, Fakira AK, et al. Identification of \n19 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nGPR83 as the receptor for the neuroendocrine peptide PEN. Sci Signal. 2016;9: ra43. \n45.  Okuda T, Haga T, Kanai Y, Endou H, Ishihara T, Katsura I. Identification and \ncharacterization of the high-affinity choline transporter. Nat Neurosci. 2000;3: 120–125. \n46.  Kiryu-Seo S, Nagata K, Saido TC, Kiyama H. New insights of a neuronal peptidase \nDINE/ECEL1: Nerve development, nerve regeneration and neurogenic pathogenesis. \nNeurochem Res. 2019;44: 1279–1288. \n47.  Philippidou P, Dasen JS. Hox genes: choreographers in neural development, architects \nof circuit organization. Neuron. 2013;80: 12–34. \n48.  Bos FL, Caunt M, Peterson-Maduro J, Planas-Paz L, Kowalski J, Karpanen T, et al. \nCCBE1 is essential for mammalian lymphatic vascular development and enhances the \nlymphangiogenic effect of vascular endothelial growth factor-C in vivo. Circ Res. \n2011;109: 486–491. \n49.  Leahy SN, Vita DJ, Broadie K. PTPN11/Corkscrew activates local presynaptic mapk \nsignaling to regulate synapsin, synaptic vesicle pools, and neurotransmission strength, \nwith a dual requirement in neurons and Glia. J Neurosci. 2024;44: e1077232024. \n50.  Ferretti R, Palumbo V, Di Savino A, Velasco S, Sbroggiò M, Sportoletti P, et al. \nMorgana/chp-1, a ROCK inhibitor involved in centrosome duplication and \ntumorigenesis. Dev Cell. 2010;18: 486–495. \n51.  Cressey R, Han MTT, Khaodee W, Xiyuan G, Qing Y. Navigating PRKCSH’s impact on \ncancer: from N-linked glycosylation to death pathway and anti-tumor immunity. Front \nOncol. 2024;14: 1378694. \n52.  Bai Y, Attardi G. The mtDNA-encoded ND6 subunit of mitochondrial NADH \ndehydrogenase is essential for the assembly of the membrane arm and the respiratory \nfunction of the enzyme. EMBO J. 1998;17: 4848–4858. \n53.  Zavodszky E, Peak-Chew S-Y, Juszkiewicz S, Narvaez AJ, Hegde RS. Identification of a \nquality-control factor that monitors failures during proteasome assembly. Science. \n2021;373: 998–1004. \n54.  Sardi SP, Murtie J, Koirala S, Patten BA, Corfas G. Presenilin-dependent ErbB4 nuclear \nsignaling regulates the timing of astrogenesis in the developing brain. Cell. 2006;127: \n185–197. \n55.  Kang J, Huang M, Li J, Zhang K, Zhu C, Liu S, et al. Enterovirus D68 VP3 targets the \ninterferon regulatory factor 7 to inhibit type I interferon response. Microbiol Spectr. \n2023;11: e0413822. \n56.  Li X, Guo H, Yang J, Liu X, Li H, Yang W, et al. Enterovirus D68 3C protease \nantagonizes type I interferon signaling by cleaving signal transducer and activator of \ntranscription 1. J Virol. 2024;98: e0199423. \n57.  Allen NJ. Astrocyte regulation of synaptic behavior. Annu Rev Cell Dev Biol. 2014;30: \n439–463. \n58.  Allen NJ, Barres BA. Neuroscience: Glia - more than just brain glue. Nature. 2009;457: \n20 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \n675–677. \n59.  Allen NJ, Eroglu C. Cell biology of astrocyte-synapse interactions. Neuron. 2017;96: \n697–708. \n60.  Chung W-S, Clarke LE, Wang GX, Stafford BK, Sher A, Chakraborty C, et al. Astrocytes \nmediate synapse elimination through MEGF10 and MERTK pathways. Nature. \n2013;504: 394–400. \n61.  Trevino AE, Müller F, Andersen J, Sundaram L, Kathiria A, Shcherbina A, et al. \nChromatin and gene-regulatory dynamics of the developing human cerebral cortex at \nsingle-cell resolution. Cell. 2021;184: 5053–5069.e23. \n62.  Choi SS, Lee HJ, Lim I, Satoh J-I, Kim SU. Human astrocytes: secretome profiles of \ncytokines and chemokines. PLoS One. 2014;9: e92325. \n63.  Sofroniew MV, Vinters HV. Astrocytes: biology and pathology. Acta Neuropathol. \n2010;119: 7–35. \n64.  Dorf ME, Berman MA, Tanabe S, Heesen M, Luo Y. Astrocytes express functional \nchemokine receptors. J Neuroimmunol. 2000;111: 109–121. \n65.  McKimmie CS, Graham GJ. Astrocytes modulate the chemokine network in a \npathogen-specific manner. Biochem Biophys Res Commun. 2010;394: 1006–1011. \n66.  Dong Y, Benveniste EN. Immune function of astrocytes. Glia. 2001;36: 180–190. \n67.  Raine CS. Morphology of myelin and myelination. Myelin. Boston, MA: Springer US; \n1984. pp. 1–50. \n68.  Simons M, Nave K-A. Oligodendrocytes: Myelination and axonal support. Cold Spring \nHarb Perspect Biol. 2015;8: a020479. \n69.  Morens DM, Folkers GK, Fauci AS. Acute flaccid myelitis: Something old and \nsomething new. MBio. 2019;10. doi:10.1128/mBio.00521-19 \n70.  Morosky S, Lennemann NJ, Coyne CB. BPIFB6 regulates secretory pathway trafficking \nand Enterovirus replication. J Virol. 2016;90: 5098–5107. \n71.  Zheng GXY, Terry JM, Belgrader P, Ryvkin P, Bent ZW, Wilson R, et al. Massively \nparallel digital transcriptional profiling of single cells. Nat Commun. 2017;8: 14049. \n72.  Hao Y, Hao S, Andersen-Nissen E, Mauck WM 3rd, Zheng S, Butler A, et al. Integrated \nanalysis of multimodal single-cell data. Cell. 2021;184: 3573–3587.e29. \n73.  Young MD, Behjati S. SoupX removes ambient RNA contamination from droplet-based \nsingle-cell RNA sequencing data. Gigascience. 2020;9: giaa151. \n \n21 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \nFigure Legends \nFigure 1. Nextstrain Phylodynamics of Enterovirus D68 Clade B. (A) Phylogenetic tree of publicly \navailable EV-D68 Clade B sequences based on VP1 region (N = 1311 sequences; collection dates: \nApril 2008 to December of 2024), highlighting the two strains used in this study - US/IL/14-18952 \n(18952) and US/MA/18-23089 (23089). (B) Temporal distribution of EV-D68 Clade B sequences, \nshown as the proportion of sequences per year.   \nFigure 2. Cell type distribution in 24-day-old hSCOs. (A) Design of scRNAseq experiment. Pooled \nhSCOs (n=12) were infected with EV-D68 strains, US/IL/14-18952 (18952) or US/MA/18-23089 \n(23089)  (10 5 PFU/group) for 48 hours at 33°C, then dissociated mechanically and enzymatically. \nCells were partitioned using the 10x NextGEM procedure and resulting libraries were sequenced \nusing short-read sequencing. Flowchart created with BioRender. (B) Uniform Manifold Approximation \nProjection (UMAP) of individual cells identi\nfied from all the samples collected in this study colored by \nassigned cell type (41,098 cells). (C) Donut plot showing the relative abundance of cell types \ncaptured across all experiments. \nFigure 3. Infection pro filing in hSCO infected with EV-D68. (A) UMAP of cellular transcriptional \nphenotypes in 24-day-old organoids infected with EV-D68 strains 18952 and 23089. Cells positive \nfor viral RNA are shown colored by viral RNA content.  (B)  Bar plots displaying the proportion of \nvRNA+ cells vs. Bystander cells in various cell types per each condition, highlighting signi ficant \ndifferences with an asterisk. (C)  Stacked bar charts showing viral read counts in each cell type. The \nnumber of total stacked cells are indicated on the top of each bar. \nFigure 4. Strain-speci fic enrichment of host pathways in infected cells. (A) Heatmap of normalized \nenrichment scores (NES) from gene set enrichment analysis (GSEA) comparing two aggregated cell \ngroups: 18952-Infected and 23089-Infected. Pathways shown represent all signi ficantly enriched \nGene Ontology Biological Process (GO-BP) terms (adjusted p < 0.05) including both positively and \nnegatively enriched gene sets. Negative NES values (blue) re flect pathways more enriched in Mock \ncells relative to the infected group. (B) Dot plots showing the expression of key genes involved in \nneuronal signaling and communication, proliferation and transcriptional activity, and suppressed \ntranslation and metabolism processes across distinct cell types in hSCOs infected with EV-D68 \nstrains 18952 and 23089. Dot size represents the percentage of cells expressing each gene, while \ncolor intensity re flects the average expression level, scaled (z-score) across all cells based on \nlog-normalized values. Midplate and VLMCs were excluded from the Mock condition, as they were \neither absent from the infected groups or detected in only one of them (e.g., VLMCs: 1 cell in \nEV-D68-18952). \nFigure 5. Strain-speci fic transcriptional responses in bystander cells. (A) Bar plot showing the \nnumber of signi ficantly di fferentially expressed genes (DEGs) between EV-D68 strains 18952 and \n23089 across bystander neurons, astrocytes, cycling astrocytes, and OPCs. DEGs were de fined as \nadjusted p-value < 0.05 and log₂FC > 0.25. B–D). Volcano plots displaying DEGs between 23089- \nand 18952-exposed bystander cells in neurons (B), astrocytes (C), and cycling astrocytes (D). Genes \nupregulated in strain 23089 are shown in green, and those upregulated in strain 18952 are in orange; \nnon-signi\nficant genes are shown in gray. Dashed lines indicate signi ficance thresholds (adjusted p = \n0.05 and log₂FC = ±0.25). (E) Bar plot showing normalized enrichment scores (NES) from GSEA of \nOPC Oligo. bystander cells, comparing transcriptional pro files between strains 23089 and 18952. \nDespite the absence of individual DEGs meeting signi ficance, GSEA revealed pathway-level \n22 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nEnterovirus D-68 Infection in Human Spinal Cord Organoids \nDábilla and Maya, et al.  \ndifferences. Bars are colored by the direction of enrichment: green for pathways enriched in 23089 \nand orange for pathways enriched in 18952. \nSupporting Information \nSupplementary Figure 1. UMAP ‘FeaturePlot’ showing Log-normalized expression of key marker \ngenes of astrocytes (SOX9 and FGFR3), cycling astrocytes (TOP2A) and neurons (MAP2 and \nTUBB3).  \nSupplementary Figure 2. Cell type characterization in hSCOs. (A) UMAP facets of all cell types in \n24-day old hSCO. (B) Donut plots showing the relative frequency of cell types in hSCOs infected with \neach EV-D68 strain and mock-infected controls.  \nSupplementary Figure 3. (A) UMAP showing the cellular transcriptional phenotypes of mGPC \nclusters in 24-days old hSCO. (B)  Heatmap of top 20 marker genes that distinguish each cluster. The \nheatmap displays scaled expression values (z-scores) of the top 20 marker genes per cluster, \ncalculated using Seurat’s ScaleData function. Values are centered and scaled per gene across all \ncells, such that 0 represents the mean expression and ±2 corresponds to approximately two \nstandard deviations above or below the mean. This highlights relative over- or underexpression \npatterns across clusters. \nSupplementary Figure 4. Scatter Plots comparing of the proportion of cells expressing ISGs \ninduced by Type I, II, and III interferons, and Volcano plots comparing the fold change in expression \nand the signi ficance, as -log 10(adjusted p-value), for comparisons of (A) Mock vs. 18952 Bystander \ncells, (B) Mock vs. 20892 Bystander cells, (C) 18952 Bystander cells vs. 20892 Bystander cells. \n23 \nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint \n\nUSC 105 and is also made available for use under a CC0 license. \n(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 \nThe copyright holder for this preprintthis version posted June 27, 2025. ; https://doi.org/10.1101/2025.06.27.661907doi: bioRxiv preprint","source_license":"Public-Domain","license_restricted":false}