Abstract The EV-A71 poses a serious threat to the health and lives of children. The EV-A71 can be transmitted by direct and indirect skin contact. Therefore, there is an urgent need to create novel skin models using human-derived cells to study the biology and pathogenesis of the virus and facilitate drug screening. Here, we use human induced pluripotent stem cells-derived skin organoids (hiPSC-SOs) as a model for EV-A71 infection and find that multiple cell types within the skin organoids, including epidermal cells, hair follicle cells, fibroblasts, and nerve cells, express EV-A71 receptors and are susceptible to EV-A71 infection. We elucidate the specific response of different cell types to EV-A71 and reveal that EV-A71 infection can degrade extracellular collagen and affect fibroblasts. We find that EV-A71 can mediate epidermal cell damage through autophagy and Integrin/Hippo-YAP/TAZ signaling pathways, thereby promoting hyperproliferation of progenitor cells. Based on this finding, we identify an autophagy-associated protein as a drug target of EV-A71 and discover an EV-A71 replication inhibitor. Altogether, these data suggest that hiPSC-SOs can be used as an infectious disease model to study skin infectious diseases, providing a valuable resource for drug screening to identify candidate virus therapeutics. Subject terms: Viral infection, Mechanisms of disease, Proteomics, Tissue engineering, Induced pluripotent stem cells __________________________________________________________________ Here, the authors establish an EV-A71 infected skin organoid model to study the virus infection pathogenesis and drug screening. The authors identify an autophagy-associated protein as drug target of EV-A71 and discover an EV-A71 replication inhibitor. Introduction Enterovirus 71 (EV-A71) is a single-stranded, positive-sense RNA virus belonging to the Enterovirus genus of the Picornaviridae family^[42]1. EV-A71 is a major pathogen of hand, foot, and mouth disease (HFMD) and mainly affects infants and young children mostly ≤ 5 years of age^[43]2. EV-A71 is primarily transmitted via the fecal-oral route. The virus enters the human body through contact with contaminated feces or hands, and the gastrointestinal tract has been identified as the primary target of infection^[44]3. The most common clinical symptoms of EV-A71 infection are low-grade fever with discomfort, macular papules or blisters on the skin of hands, feet, and buttocks, and painful ulcerative lesions in the throat, mouth, and tongue^[45]4. Although the virus is easily cultured from skin vesicles and oral secretions, the replication site and pathogenesis of EV-A71 in human skin remain unknown. In addition, EV-A71 is a highly neurotropic virus, and the central nervous system (CNS) is also an important target for EV-A71 infection. Following EV-A71 infection of the CNS, inflammation can occur in different anatomical regions or simultaneously in multiple regions, causing encephalitis, meningitis, neurogenic pulmonary edema, and acute flaccid paralysis^[46]5. Thus, timely identification and treatment of EV-A71 infections is important. In the past decade, EV-A71 has emerged as a serious threat to public health in the Asia–Pacific region^[47]6–[48]10. Unfortunately, currently, there are no specific antivirals available to clinically treat this infection. Appropriate cell and animal models are needed to better understand EV-A71 pathogenesis and to aid the development of effective drugs. Until now, none of the EV-A71 animal models have been able to recapitulate all aspects of the human disease. These models may show diversity in the tissue tropisms, primary viral replication sites, and disease manifestations of EV-A71 infections with human^[49]11. Therefore, there is an urgent need to develop new techniques to discover new targets and evaluate potential therapeutics. Human pluripotent stem cells have the ability to self-organize into complex organoids in vitro, and these 3D organoids mimic the structure and function of natural human organs and are capable of modeling in vivo organ development and disease^[50]12. In recent years, multiple human organoid models have been generated to establish novel 3D platforms for in vitro disease research, especially SARS-CoV-2-infected organoid models^[51]13–[52]16. In vitro models are crucial for understanding disease pathogenesis and developing therapeutic targets to mitigate the harmful impact of infectious diseases on public health. Since traditional 2D cellular models cannot accurately simulate complex microenvironments at the organ level and in vivo animal models exhibit species differences, human stem cell-derived organoids are excellent alternatives^[53]17. Skin organoids generated from human pluripotent stem cells are composed of stratified epidermis, fat-rich dermis, pigmented hair follicles, sebaceous glands, and nerve bundles consisting of neurons and Schwann cells^[54]18. Skin organoids have been applied in microbial infection modeling of SARS-CoV-2, Mpox virus^[55]19 and Staphylococcus aureus, showing their strong potential as a tool to investigate the mechanisms of viral infection and assess drug efficacy^[56]14,[57]20. Therefore, using skin organoids can be a powerful strategy for studying the effects of EV-A71 infection on the skin and nervous system. In this study, human induced pluripotent stem cell-derived skin organoids (hiPSC-SOs) containing epithelial cells, hair follicles and nerve cells were established to study the effects of EV-A71 infection on the skin. Single-cell RNA sequencing (scRNA-seq) and proteomics profiling of EV-A71-infected skin organoids were first constructed to study the virulence, molecular characteristics and pathogenesis of EV-A71 infection. We identified fibroblasts and epidermal cells as the main target cells of EV-A71 and elucidated their roles in the abnormal pathological phenotypes caused by viral infections. Autophagy and aberrant proliferation of epidermal progenitor cells (PCs) were found to be important mechanisms for viral infection of epidermal cells, which may be associated with tissue tumorigenesis. By further analysis of proteins related to viral replication, a potential drug target was discovered. Results hiPSC-SOs are susceptible to EV-A71 infection We differentiated hiPSCs into skin organoids using an improved method according to the previously reported stepwise protocol^[58]14 (Supplementary Fig. [59]1A). Immunofluorescence staining analysis results validated the expression of the basal stem cell marker keratin 5 (KRT5), hair follicle stem cell markers keratin 15 (KRT15) and keratin 17 (KRT17), as well as the pan-neuronal marker tubulin β class III (TUJ1), in the hiPSC-SOs (Supplementary Fig. [60]1B). Electron microscopy showed the different hierarchical structures of hair follicles, secreted melanosomes, extracellular matrix, and desmosomes in the epidermis (Supplementary Fig. [61]1C). Next, to study the different cell types that are affected by EV-A71 infection, the hiPSC-SOs were subjected to EV-A71 infection and the proteomics analysis was performed on the infected organoids (Fig. [62]1A). With increasing EV-A71 infection time, the amount of viral replication in the hiPSC-SOs gradually increased, indicating that the hiPSC-SOs were susceptible to EV-A71 infection (Fig. [63]1B, C and Supplementary Fig. [64]2A). Whole-mount immunostaining of EV-A71-infected hiPSC-SOs revealed large amounts of virus protein appeared at different areas of hair follicles (KRT17, KRT82, and KRT71) and nerves (TUJ1 and S100 calcium binding protein B (S100B)) in organoids (Fig. [65]1D and Supplementary Fig. [66]2B), suggesting that the virus may infect different types of cells in hiPSC-SOs. Immunostaining analysis of organoid sections further validated that EV-A71 protein co-localized with the markers of hair follicle cells (KRT17, KRT82, and KRT71), as well as neurons (TUJ1) and S100B+ Schwann-like cells (Fig. [67]1E and Supplementary Figs. [68]2C, [69]S3), indicating the virus may infect these types of cells. Fig. 1. hiPSC-SOs are susceptible to EV-A71 infection. [70]Fig. 1 [71]Open in a new tab A Workflow of EV-A71 virus-infected human skin organoids (hiPSC-SOs) and proteomics analysis. B Immunofluorescence of EV-A71 virus protein in the infected hiPSC-SOs (scale bar: 50 μm). C Evaluation of the EV-A71 viral infection efficiency at 2, 4, and 8 days post-infection (dpi). The quantification analysis is performed from three images for each group. The data are shown as the mean ± SD (n = 3 per group). Significant differences are determined by unpaired two-tailed t test. D Whole-mount staining of EV-A71 virus protein KRT17, and TUJ1 in one side and the other side of infected hiPSC-SOs at 8 dpi (scale bar: 200 μm). E Immunofluorescence of EV-A71 virus protein with different cellular markers, KRT71, KRT82, TUJ1, and S100, in infected hiPSC-SOs at 8 dpi (scale bar: 20 μm). The white triangles indicate where the viruses co-localize with corresponding cellular markers. The immunofluorescence experiments are repeated three times. hiPSC-SOs: human induced pluripotent stem cell-derived skin organoids. Source data are provided as a Source Data file. Then, the quantitative proteomics analyses were performed on EV-A71-infected hiPSC-SOs over time and a total of 3128, 3198, 3124, and 3053 proteins were identified in hiPSC-SOs at 0 (control), 2, 4 and 8 days post-infection (dpi) (Supplementary Data [72]1). The proteome profiling showed that the protein expression pattern of the hiPSC-SOs without EV-A71 infection was quite different from that of the hiPSC-SOs infected with EV-A71 (Supplementary Fig. [73]4A, B). The pairwise comparisons were performed by a moderated t test to identify the differentially expressed proteins (DEPs) between the hiPSC-SOs with EV-A71 infection at 2 (n = 3), 4 (n = 3), and 8 (n = 3) dpi and the hiPSC-SOs without EV-A71 infection (control group, n = 3). The DEPs between EV-A71-infected and control samples were determined based on a Benjamin–Hochberg (BH) adjusted p value < 0.01 and log2 (EV-A71 infected/Control) > 0.585 (upregulated) or < -0.585 (downregulated). Among these DEPs, we found 552, 335, and 371 upregulated proteins and 187, 244, and 364 downregulated proteins in the hiPSC-SOs infected with EV-A71 for 2, 4, and 8 days, respectively, compared to the hiPSC-SOs without EV-A71 infection (Fig. [74]2A, Supplementary Data [75]2). Biological process enrichment analysis of all the DEPs showed that the upregulated proteins were enriched in mRNA processing, translational initiation, viral entry into host cells, and DNA replication (Supplementary Fig. [76]4C), indicating that virus cycle-associated processes may be activated in EV-A71-infected hiPSC-SOs. The downregulated proteins were enriched in collagen fibril organization, basement membrane organization, keratinization, and fibroblast proliferation (Supplementary Fig. [77]4C), indicating that the abundance of skin cell- and extracellular matrix-related proteins was severely downregulated in EV-A71-infected hiPSC-SOs. Fig. 2. Quantitative expression profiles and biological characteristics of proteins identified in EV-A71-infected hiPSC-SOs. [78]Fig. 2 [79]Open in a new tab A Volcano plots of –log10 p value vs. log2 protein abundance comparisons between the EV-A71-infected hiPSC-SOs (2, 4, and 8 dpi) and control hiPSC-SOs (0 dpi) (n = 3 for each group). The pairwise comparisons are performed by the two-sided moderated t test to identify the DEPs between the hiPSC-SOs with EV-A71 infection at 2 (n = 3), 4 (n = 3), and 8 (n = 3) dpi and the hiPSC-SOs without EV-A71 infection (0 dpi, n = 3). The DEPs between EV-A71-infected and control samples are determined based on a BH adjusted p value < 0.01 and log2 (EV-A71 infected/Control) > 0.585 (upregulated) or < -0.585 (downregulated), and labeled in red (upregulated) or blue (downregulated), respectively. B Biological process enrichment analysis of proteins that are specifically highly and lowly expressed in EV-A71-infected hiPSC-SOs at 2 (n = 3), 4 (n = 3), and 8 (n = 3) dpi compared to the controls (0 dpi, n = 3). The specific highly and lowly expressed proteins are determined with a BH -adjusted p value < 0.01 and log2FC > 0.585 or < -0.585, respectively. The fold change (FC) represents the ratio of the normalized intensity of the protein identified in the hiPSC-SOs at the indicated infection time to those at all other samples. The enrichment analysis is performed using the one-sided hypergeometric test and the p values are provided in source data. The enriched terms of high and low expressed proteins are indicated in red and blue, respectively. The size of the dot indicates the number of proteins belonging to each term. The color scale indicates the enrichment p value. C Expression heatmap of the proteins that are involved in specific biological functions of EV-A71-infected hiPSC-SOs at 0, 2, 4, and 8 dpi. The red and blue boxes indicate proteins with increased and decreased abundance levels, respectively. The clusters of proteins associated with similar biological processes are grouped and labeled. BH: Benjamin–Hochberg; DEP: differentially expressed protein; dpi: days post-infection; hiPSC-SOs: human induced pluripotent stem cell-derived skin organoids. Source data are provided as a Source Data file. The above results suggest that hiPSC-SOs can be used as in vitro models of EV-A71 infection. Time-dependent molecular features of EV-A71-infected hiPSC-SOs To investigate the dynamic changes in EV-A71-infected hiPSC-SOs over time, we examined the proteomics profiling at 2, 4 and 8 dpi (Supplementary Fig. [80]4A, B). The specifically highly and lowly expressed proteins at different infection timepoint were defined as those with BH adjusted one-way ANOVA p-value < 0.01 and log2FC > 0.585 or < -0.585, respectively. The fold change (FC) represents the ratio of the normalized intensity of the protein identified in the hiPSC-SOs at the indicated infection time to those in all other samples. The results showed specific high and low expression patterns of protein abundances at the assessed time points (Supplementary Fig. [81]4D, Supplementary Data [82]3). Four high expression patterns of protein abundances were found at 0, 2, 4, and 8 days after EV-A71 infection. The functional annotation of these proteins revealed the proteomic signatures at each timepoint. At 2 dpi (early stages of infection), the proteins with high abundance were mainly functionally enriched in virus processes, such as response to virus and viral genome replication, as well as DNA replication, such as prereplicative complex assembly involved in nuclear cell cycle DNA replication and DNA replication initiation, indicating that in the early stages of hiPSC-SO EV-A71 infection, the main biological processes of host cells appear to help virus replication (Fig. [83]2B). For example, the DEAD box protein DDX21, an RNA helicase, promotes ribosomal RNA processing and transcription from polymerase II 1 and 2 and is highly expressed at 2 dpi (Fig. [84]2C). Glycolytic processes (ALDOA, ALDOC, ENO2, GPI, and HK2) and glucose import across the plasma membrane (SLC2A1 and SLC2A3) seem to occur at 4 dpi (Fig. [85]2B, C), thereby suggesting enormous cellular energy expenditure at this point during the infection. Finally, the proteins with high abundance were enriched in immune responses, such as cell chemotaxis (CCL20, CXCL12, DEFB103A, and TPBG), responses to cytokines (ACP5, PTGES, and TIMP1), and responses to hypoxia (ATP1B1, ENG, and PLOD2), at 8 dpi, thereby indicating that immune and hypoxic responses occur in the later stages of EV-A71 infection (Fig. [86]2B, C). Similarly, we obtained the specific low expression patterns of protein abundances via expression profile clustering (Supplementary Fig. [87]4D). Hair follicle morphogenesis (KRT25, KRT27, KRT34, and KRT71) and hair differentiation (KRT33A, KRT36, and KRT40)-associated proteins seemed to be downregulated at 4 dpi (Fig. [88]2B, C), indicating that hair follicle may be an important structure and its function would be affected after EV-A71 infection. Immunostaining analysis further validated that EV-A71 can directly infect hair follicle cells at the different timepoint assessed (Supplementary Fig. [89]4E). We further observed that proteins involved in epidermal development, such as basal stem cell markers (KRT5 and KRT15) and epithelium markers (IVL and KRT12) were significantly dysregulated at 8 dpi in EV-A71-infected hiPSC-SOs (Fig. [90]2B, C). In addition, the main components of the basement membrane (COL7A1, COL17A1, LAMA1, LAMB1, LAMC1, NID1, and NID2), the basal stem cell markers, the specialized anchoring complexes of skin, hemidesmosomes (ITGB4 and ITGA6) and adhesion receptors (ITGA3) were also downregulated significantly at 8 dpi in EV-A71-infected hiPSC-SOs (Fig. [91]2B, C). The ECM organization-associated proteins such as collagen fibril-related proteins (COL1A1, COL1A2, COL2A1, COL3A1, COL5A1, COL11A1, COL14A1, and LUM) and elastic fiber assembly-related proteins (LOX and EMILIN1) were also downregulated in EV-A71-infected hiPSC-SOs at 8 dpi (Fig. [92]2B, C). In conclusion, the above data showed the processes of dynamic interactions of EV-A71 with hiPSC-SOs and indicated that the viral infection can lead to the disorder of specific physiological structures and functions of the skin. EV-A71 infection depletes reticular fibroblasts Next, we performed scRNA-seq of EV-A71-infected hiPSC-SOs to investigate the cell tropism of EV-A71 and its impact on the different cell types within hiPSC-SOs (Fig. [93]3A, Supplementary Data [94]4). scRNA-seq data of the hiPSC-SOs revealed the presence of fibroblast (COL5A1 and COL6A3), nerve cell (HES1 and PTN), keratinocyte (KRT5 and KRT14), melanocyte (PMEL and TYRP1), merkel cell (KRT8 and TTR), chondrocyte (CHI3L2 and SOX5), and cycling progenitor cell (TOP2A and CDK1) (Fig. [95]3B, Supplementary Data [96]4). The key factors involved in EV-A71 entry into host cells, including ANXA2^[97]21, FN1^[98]22, NCL^[99]23, PHB2^[100]24, PPIA^[101]25, SCARB2^[102]26, and VIM^[103]27, were found to be expressed in main types of cells identified in scRNA-seq data (Fig. [104]3C). Correspondingly, the numbers of fibroblasts, nerve cells, and epidermal cells showed greatest changes before and after EV-A71 infection (Fig. [105]3D). Fig. 3. EV-A71 infection depletes reticular fibroblasts. [106]Fig. 3 [107]Open in a new tab A UMAP plots of scRNA-seq data generated from EV-A71-infected hiPSC-SOs at 0 and 8 dpi. A total of 33,971 cells are represented. The major cell groups are manually annotated and labeled with different colors. FB: fibroblast; NC: nerve cell; EpC: epidermal cell; Mel: melanocyte; Mer: merkel cell; Cho: chondrocyte; PC: progenitor cell. B Expression and percentage distribution of the key gene markers of different cell clusters. The statistical significance of cell type marker genes is determined using the two-sided Wilcoxon rank-sum test with Bonferroni-adjusted p value < 0.05. C Gene expression levels of key receptors (ANXA2, FN1, NCL, PHB2, PPIA, SCARB2, and VIM) involved in EV-A71 entry in scRNA-seq data of FB (n = 22,649), NC (n = 3,712), EpC (n = 4,779), Mel (n = 1,077), Mer (n = 639), Cho (n = 336), and PC (n = 599). Boxes represent the median values and the first and third quartiles. The upper and lower whiskers represent ranges that extending up to 1.5 the interquartile range. The outliers with the values below the first quartile or above the third quartile are represented with stars. D The proportion of each cell type in normal control and EV-A71-infected hiPSC-SOs. E Biological process enrichment analysis of upregulated genes between EV-A71-infected and normal hiPSC-SOs in proliferative fibroblast and reticular fibroblast. The DEGs between EV-A71-infected and control hiPSC-SOs are determined based on a BH-adjusted p value of <0.05 and log2 (EV-A71 infected/Control) > 0.25 (upregulated). The enrichment analysis is performed using the one-sided hypergeometric test and the p values are provided in source data. F Mechanism diagram of EV-A71 infecting reticular fibroblasts. When reticular fibroblasts are infected by viruses, TNF-mediated signaling pathway (TNFR1, CREB3, and C/EBPβ) is activated, producing pro-inflammatory factors (IL6 and CXCL8) and MMPs (MMP1, MMP3, and MMP13) that metabolize extracellular collagen. Ultimately, skin inflammation and collagen loss combine to cause the onset of skin aging. BH, Benjamin–Hochberg; DEG: differentially expressed gene; dpi: days post-infection; hiPSC-SOs: human induced pluripotent stem cell-derived skin organoids. Source data are provided as a Source Data file. As the largest population in target cells of EV-A71, we found that the number of fibroblasts decreased after virus infection in hiPSC-SOs (Fig. [108]3D and Supplementary Fig. [109]5A). Next, we found that fibroblasts consisted of different cell types, such as proliferative fibroblast (FGF7 and NRN1), reticular fibroblast (MFAP5 and COL11A1), MKI67 PC (CENPW and F2R), and chondrocyte (COL2A1 and ACAN) (Supplementary Fig. [110]5B–D, Supplementary Data [111]4). Among them, the number of reticular fibroblasts, which are a type of fibroblast that synthesizes collagen alpha-1 (III) to produce and maintain the thin fibrous networks that are the foundation of most lymphoid organs, markedly decreased after EV-A71 infection (Supplementary Fig. [112]5E, F), indicating that dermal tissue lost the ability to provide structural support after EV-A71 infection and that the balance of the immune response was disrupted. The biological process enrichment analysis was performed on the upregulated differentially expressed genes (DEGs) in reticular fibroblasts and proliferative fibroblasts between the EV-A71-infected hiPSC-SOs (8 dpi) and control hiPSC-SOs (0 dpi). The results showed that both reticular fibroblasts and proliferative fibroblasts exhibited a significant enrichment in host-virus interaction, DNA repair, apoptosis, and aging related proteins (Fig. [113]3E). Reticular fibroblasts were enriched in more processes related to host of viral transcription and response to virus, indicating that reticular fibroblasts are the main target of EV-A71 virus infection (Fig. [114]3E). In addition, immune response pathways, such as TNF-α and IL-1 mediated signaling pathways, were enriched in EV-A71-infected reticular fibroblasts (Fig. [115]3E), leading to the upregulation of IL-6 and CXCL8 (Fig. [116]3F). Activation of the TNF pathway (TNFR1, CREB3, and C/EBPβ) leads to upregulation of matrix metalloproteinases (MMPs), which can degrade several collagens^[117]28,[118]29. Results showed that MMP1, MMP3, and MMP13 were upregulated in the EV-A71-infected reticular fibroblasts (Fig. [119]3F), which may be the key pathway for the massive reduction of extracellular matrix components such as collagen in previous proteomic results (Figs. [120]2C, [121]3F). Loss of dermal collagen accompanied by the production of pro-inflammatory factors will contribute to skin senescence^[122]30. These results suggested that EV-A71 infection depletes fibroblast subsets and spreads senescence and thus may trigger skin aging. EV-A71 infection promotes NNMT + SFRP1+ progenitor cells’ proliferation through the Hippo-YAP/TAZ signaling pathway To further investigate the impact of EV-A71 on epidermal cells, we analyzed the scRNA-seq data of epidermal cell subtypes including basal stem cell, mature epidermal cell (mature EpC), and NNMT + SFRP1+ progenitor cell (NNMT + SFRP1 + PC) (Fig. [123]4A–D, Supplementary Data [124]4). Results showed that genes associated with host-virus interaction, mRNA processing, positive regulation of transcription, viral translational termination-reinitiation, aging, and the innate immune response were upregulated in these three epidermal cell types (Fig. [125]4E), indicating disruption of epidermal cell homeostasis by viral infection. In addition, genes associated with skin development, skin barrier establishment, collagen fibril organization, elastic fiber assembly, and tissue regeneration were downregulated in epidermal cells (Fig. [126]4E), which is consistent with our proteome results, indicating that EV-A71 infection of hiPSC-SOs affects the growth and barrier function of epidermal cells, as well as their matrix environment in the dermis. Interestingly, we found the number of NNMT + SFRP1+PCs nearly doubled after virus infection (Fig. [127]4D). Immunofluorescence staining confirmed the surge of NNMT + SFRP1+PCs in a time-dependent manner after viral infection (Fig. [128]4F). Several upregulated genes associated with viral processes, such as antiviral defense, viral genome replication, defense response to virus were enriched in the NNMT + SFRP1+PCs of the EV-A71-infected hiPSC-SOs, compared to the control group (Fig. [129]4E). KEGG pathway analysis showed that the upregulated genes enriched in the NNMT + SFRP1+PCs were associated with EGF/EGFR signaling (MAP3K2, MAP2K2, etc.), TGF-beta signaling (MMP1, TNC, etc.), integrin-mediated cell adhesion (ITGB1, ITGAV, etc.), tight junction (ACTR2, ACTR3, etc.), as well as Hippo signaling (FZD8, TAZ, etc.) pathways (Fig. [130]4G). Activation of the Integrin and Hippo signaling pathways in NNMT + SFRP1+PCs promotes YAP/TAZ expression, which may be an important pathway leading to increased NNMT + SFRP1 + PC numbers after EV-A71 infection (Fig. [131]4H). The PC marker NNMT plays a key role in the development and progression of several cancers^[132]31. Therefore, in addition to the aberrant proliferation of PCs, activation of YAP/TAZ could also induce cancer cell phenotypes of PCs^[133]32. The above results suggested that EV-A71 infection may promote NNMT + SFRP1+PCs’ proliferation through the Hippo-YAP/TAZ signaling pathway. Fig. 4. EV-A71 infection promotes NNMT + SFRP1 + PCs’ proliferation through the Hippo-YAP/TAZ signaling pathway. [134]Fig. 4 [135]Open in a new tab A UMAP plots showing four subpopulations of epidermal cells in EV-A71-infected hiPSC-SOs at 0 and 8 dpi. B Dot plot showing the expression of representative marker genes for NNMT + SFRP1 + PC, basal stem cell, and mature EpC. The statistical significance of cell type marker genes is determined using the two-sided Wilcoxon rank-sum test with Bonferroni-adjusted p value of <0.05. C Expression distribution of representative marker genes for NNMT + SFRP1 + PC (NNMT and SFRP1), basal stem cell (TP63 and COL17A1), and mature EpC (KRT1 and IVL). D Proportions of different types of epidermal cell in EV-A71-infected and normal hiPSC-SOs. E Biological process enrichment analysis of the DEGs in NNMT + SFRP1 + PC, basal stem cell, and mature EpC between EV-A71-infected and normal hiPSC-SOs. The enriched terms of upregulated and downregulated DEGs are indicated in red and blue, respectively. The size of the dot indicates the number of genes belonging to each term. The color scale indicates the enrichment p value. F Immunofluorescence of the EV-A71 virus protein and NNMT + SFRP1 + PC markers (NNMT and SFRP1) in the EV-A71-infected hiPSC-SOs at 0, 4, and 8 dpi (scale bar: 100 μm). The bar charts quantify the expression levels of NNMT + SFRP1+PCs’ markers NNMT and SFRP1 in the EV-A71-infected hiPSC-SOs with infection time. The data are shown as the mean ± SD (n = 3 per group). G KEGG pathway analysis of the upregulated genes in the NNMT + SFRP1 + PC of EV-A71-infected hiPSC-SOs at 8 dpi compared with normal controls. The enrichment analysis is performed using the one-sided hypergeometric test and the p values are provided in source data. H Mechanism diagram of EV-A71 infection-induced proliferation of NNMT + SFRP1 + PC. When NNMT + SFRP1+PCs are infected by viruses, intergrin-mediated signaling pathway (ITGAV and ITGB1) and following Hippo signaling pathway (YAP, TAZ, TEAD1, and BIRC2) are activated to promote proliferation and invasion of NNMT + SFRP1+PCs. The immunofluorescence experiments are repeated three times. dpi: days post-infection; EpC: epidermal cell; hiPSC-SOs: human induced pluripotent stem cell-derived skin organoids; PC: progenitor cell. Source data are provided as a Source Data file. Collectively, these results suggested that EV-A71 mediates abnormal epidermal phenotypes by simultaneously disrupting the normal ecological balance of epidermal cells and inducing overproliferation of PCs, which may be associated with tumorigenesis. Drug blocks EV-A71 replication Next, in order to deeply investigate the reasons for the increased proliferation of NNMT + SFRP1+PCs after EV-A71 infection and screen candidate antiviral targets, we performed the Wilcoxon rank-sum test to identify DEGs in NNMT + SFRP1+PCs from the scRNA-seq data between the EV-A71-infected hiPSC-SOs and control samples. The genes with Bonferroni-adjusted p value < 0.05 and log2 (EV-A71 infected/Control) > 0.25 or < -0.25 are considered as upregulated or downregulated DEGs. There were 744 upregulated and 472 downregulated DEGs in the hiPSC-SOs infected with EV-A71 compared to those without EV-A71 infection (Supplementary Fig. [136]6A). Within the upregulated DEGs in NNMT + SFRP1+PCs, 12 genes were overlapped with the proteins that were simultaneously upregulated in EV-A71-infected hiPSC-SOs at 2, 4, and 8 dpi (Supplementary Fig. [137]6B and Fig. [138]5A). Among them, HMGB1 (high mobility group box 1), a protein associated with autophagy^[139]33 and viral infection^[140]34,[141]35, aroused our interest. Furthermore, we detected the gene expression level of HMGB1 based on scRNA-seq data and found that it mainly expressed in fibroblasts, nerve cells, epidermal cells, and hair follicle cells, which was consistent with the expression levels of some of the EV-A71 receptor (Supplementary Fig. [142]6C and Fig. [143]3C). In addition, HMGB1 showed increased expression levels in all cell types after EV-A71 infection; then the knockdown experiments were performed on HMGB1. The results showed knockdown of HMGB1 can prevent EV-A71 replication in cell lines (Fig. [144]5B and Supplementary Fig. [145]6D). Fig. 5. NSC167409 blocks EV-A71 replication both in vivo and in vitro. [146]Fig. 5 [147]Open in a new tab A Overlap of the upregulated DEGs in NNMT + SFRP1+PCs with the upregulated DEPs in the EV-A71-infected hiPSC-SOs at 2, 4, and 8 dpi compared to 0 dpi. The red boxes present the log2 average expression of genes in NNMT + SFRP1+PCs between the EV-A71-infected hiPSC-SOs and controls. B Western blot analysis of HMGB1, the EV-A71 virus protein (EV-A71-VP1), and tubulin in the siHMGB1-1, siHMGB1-2, siControl, and normal control (NC) samples of Vero cells infected with EV-A71 (1 MOI) for 48 h. The quantitative data are shown as the mean ± SD (n = 3 per group). Significant differences are determined by unpaired two-tailed t test and the p values are provided in source data. C Western blot analysis of the EV-A71 virus protein in EV-A71-infected (1 MOI) Vero cells that treated with NSC167409. The samples are treated with NSC167409 immediately after being infected with EV-A71 and collected at 30 hpi for analysis. D Antiviral activity and cytotoxicity analysis of NSC167409 on RD and Vero cells. The samples are treated with NSC167409 immediately after being infected with EV-A71 (100 TCID[50]/ml) and collected at 72 hpi for analysis. The antiviral activity and cytotoxicity of NSC167409 are measured using CellTiter-Glo cell viability assay kit. The data are shown as the mean ± SD (n = 3 per group). E Immunofluorescence of EV-A71 virus protein and KRT5 in the EV-A71-infected hiPSC-SOs at 0 and 8 dpi with or without NSC167409-treatment (scale bar: 500 μm). F Workflow of in vivo transplantation of hPSC-SOs, drug administration, and EV-A71 infection on mice. The prophylactic treatment is performed on mice with NSC167409 (2 mg/kg) by subcutaneous injection four hours before EV-A71 infection. G Immunofluorescence of the EV-A71 virus protein and KRT5 in the EV-A71-infected hiPSC-SO xenografts at 0 and 8 dpi with or without NSC167409-treatment (scale bar: 100 μm). The immunofluorescence experiments are repeated three times. DEG: differentially expressed gene; DEP: differentially expressed protein; dpi: days post-infection; hiPSC-SOs: human induced pluripotent stem cell-derived skin organoids; hpi: hours post-infection. Source data are provided as a Source Data file. Moreover, we found that the expression level of HMGB1 increased in a time-dependent manner with the increase of virus infection time (Supplementary Fig. [148]7A), consistent with previous proteomic and scRNA-seq results. To identify drug candidates capable of inhibiting EV-A71 replication, cell lines and hiPSC-SOs were treated with an FDA-approved HMGB1 inhibitor (NSC167409, 0.3-10 μM). The cells were treated with NSC167409 both immediately (Fig. [149]5C, D) and several hours (Supplementary Fig. [150]7B) after being inoculated with EV-A71; the results showed that the inhibition rate of EV-A71 levels increased both in rhabdomyosarcoma (RD) cells and Vero cells with a dose-dependent manner of drug treatment (Fig. [151]5C, D and Supplementary Fig. [152]7B). The inhibitory effect of NSC167409 on EV-A71 was also observed in hiPSC-SOs; and the drug showed a more significant inhibitory effect when treated at 24 hours post-infection (hpi) than at 72 hpi (Supplementary Fig. [153]7C). The immunofluorescence results further validated the inhibitory effect of NSC167409 on EV-A71 replication (Fig. [154]5E). To evaluate the drug activities in vivo, we performed the prophylactic treatment on humanized mice carrying hiPSC-SOs with NSC167409 (2 mg/kg) before intraxenograft inoculation with the EV-A71 virus (Fig. [155]5F). EV-A71 staining was significantly decreased in the epidermal cells of hiPSC-SO xenografts from humanized mice treated with NSC167409 and EV-A71, suggesting that the drug successfully prevented the virus from multiplying in vivo (Fig. [156]5G). The above results indicated that NSC167409 can effectively inhibit virus replication both in vitro and in vivo. EV-A71 mediates epidermal cell dysfunction via autophagy pathways in hiPSC-SOs To further investigate the specific role of HMGB1 in EV-A71 infection, we analyzed HMGB1 expression levels in different cellular subpopulations before and after EV-A71 infection. We found that HMGB1 was highly expressed in the NNMT + SFRP1+PCs, mature EpCs, basal stem cells, proliferative fibroblasts, and reticular fibroblasts in scRNA-seq data of EV-A71-infected hiPSC-SOs compared to the control group (Fig. [157]6A), indicating that HMGB1 may be related to EV-A71 infection in these cells. Among them, we found that HMGB1 expression level was highest in NNMT + SFRP1+PCs and proliferative fibroblasts (Fig. [158]6A). Since the number of NNMT + SFRP1+PCs increased specifically after viral infection (Figs. [159]4D, [160]6A), suppression of HMGB1 with NSC reduced the number of NNMT + SFRP1+PCs, indicating that the specific high expression level of HMGB1 was closely related to the number of NNMT + SFRP1+PCs (Supplementary Fig. [161]7D). Fig. 6. EV-A71 mediates epidermal cell dysfunction via autophagy pathways in hiPSC-SOs. [162]Fig. 6 [163]Open in a new tab A Violin plots illustrating the expression of HMGB1 in selected subclusters of epidermal cell and fibroblast of the EV-A71-infected hiPSC-SOs at 8 dpi. B Pathway enrichment analysis of DEGs in NNMT + SFRP1 + PC, basal stem cell, and mature EpC in the EV-A71-infected hiPSC-SOs between 8 dpi and 0 dpi. The upregulated and downregulated DEGs are indicated in red and blue, respectively. The size of the dot indicates the percentage of cells with gene expression in a cell cluster, while the color indicates the average expression level of the gene. The DEGs of the EV-A71-infected hiPSC-SOs between 8 dpi and 0 dpi are determined based on the two-sided Wilcoxon rank-sum test with Bonferroni-adjusted p value of <0.05 and log2 (EV-A71 infected/Control) > 0.25 (upregulated) or < −0.25 (downregulated). C Immunofluorescence of the EV-A71 virus protein, NNMT, TAZ, YAP, TEAD1, LC3B, and HMGB1 in the EV-A71-infected hiPSC-SOs at 0, 4 and 8 dpi with or without NSC167409-treatment (scale bar: 10 μm). The experiments are repeated three times. D Mechanistic diagram of the role of HMGB1 in EV-A71 infection of NNMT + SFRP1+PCs. HMGB1 positively correlates with viral replication. HMGB1 promotes autophagy (LC3B) and the production of pro-inflammatory factors (IL6, IL8, and MMP1), and enhances YAP/TAZ signaling to promote the NNMT + SFRP1+PCs’ proliferation. When an inhibitor (NSC167409) is used to inhibit the expression of HMGB1, not only is viral replication and its concomitant inflammation suppressed, but autophagy and YAP/TAZ signaling are also suppressed. DEG: differentially expressed gene; dpi: days post-infection; hiPSC-SOs: human induced pluripotent stem cell-derived skin organoids; EpC: epidermal cell; FB: fibroblast; PC: progenitor cell. Source data are provided as a Source Data file. By gene enrichment analysis, we found that several signaling pathways were activated, such as HIF-1 (EGLN3, EIF4BP1, ALDOC, ALDOA, EGLN1, CUL2, SLC2A1, EGFR, VEGFA, and ELOB) and autophagy (HMGB1, HIF1A, TANK, VMP1, RB1CC1, DDIT4, and NRBF2) in epidermal cells in EV-A71-infected hiPSC-SOs (Fig. [164]6B), indicating that EV-A71 replication and EV-A71-host interactions may be achieved through HIF-1 and autophagy signaling pathways. Immunostaining analysis further validated that over time, EV-A71 infection can upregulate the expression of microtubule-associated protein 1A/1B light chain 3B (LC3B), which is associated with the formation of autophagosomes, accompanied by an increase of the viral load in mature EpCs, basal stem cells and NNMT + SFRP1+PCs of the EV-A71-infected iPSC-SOs, as well as the humanized mice carrying hiPSC-SOs with EV-A71 infection (Fig. [165]6C and Supplementary Fig. [166]7A, E, F, [167]8A). Moreover, the expression of the autophagy inhibitor BCL2 was downregulated, further validating the omics analysis results. Suppression of HMGB1 with NSC reduced the expression of LC3B in vitro and vivo (Fig. [168]6C and Supplementary Fig. [169]8B), indicating that EV-A71 may replicate through the autophagy pathway. In previous study, we found EV-A71 infection promotes NNMT + SFRP1+PCs’ proliferation through the Hippo-YAP/TAZ signaling pathway. Results also showed that YAP, TAZ, and TEAD1 of Hippo-YAP/TAZ signals were highly expressed in a time-dependent, accompanied by an increase of the viral load in EV-A71-infected organoids (Fig. [170]6C). Further, suppression of HMGB1 with NSC can reduced the expression of them, accompanied by decreased expression of LC3B and HMGB1 and the number of NMT + SFRP1+PCs (Fig. [171]6C). Following inhibition of YAP/TAZ signaling by 5 μΜ Verteporfin, we found that the expression of NNMT was significantly decreased in infected organoids (Supplementary Fig. [172]9). These results indicated that the increase in the number of NNMT + SFRP1+PCs may be related to autophagy-YAP signaling (Fig. [173]6D), which is also the signaling pathways in cancer progression^[174]36,[175]37. Gene enrichment analysis also showed that immune response-related pathways, such as the IL-6 signaling pathway (TIMP1, PRDM1, NR2F6, IL6, and IL6ST), IL-2 signaling pathway (MAP2K2, SHC1, RPS6KB2, PTPN11, CBL, HRAS, FYN, SOS1, and JAK1), IL-17 signaling pathway (CXCL8, MMP1, MMP3, MMP13, S100A9, S100A8, S100A7, S100A7A, JUND, MAPK6, and HSP90B1) and TGFβ signaling pathway (MAP2K3, ROCK1, ITGA2, KLF6, ZFYVE16, COPS5, SNW1, SKIL, PAK2, ITGB1, TNC, DAB2, and ZEB1), were activated in EV-A71-infected hiPSC-SOs (Fig. [176]6B), indicating that EV-A71 infection may induce inflammation in epidermal cells. Immunofluorescen indicated that HMGB1 inhibitor also attenuated inflammatory response caused by viral infection due to the downregulated expression levels of IL6, IL8, MMP1, and MMP3 in the EV-A71-infected xenografts with NSC treatment in mice, compared to those without NSC treatment (Supplementary Fig. [177]10). These results suggested that by suppressing the HMGB1, epidermal cell dysfunction and inflammatory responses caused by viral infections were effectively alleviated. Our data illustrated an important mechanism for NSC as a potential drug for inhibiting epidermal PC abnormal proliferation by anti EV-A71 replication. Discussion The human stem cell-derived organoid model is an excellent in vitro model for studying viral tropism and pathogenesis as it mimics the physiological organ microenvironment. Previous studies have used organoid models to study the infection characteristics of Zika^[178]38 and SARS-CoV-2^[179]13 and potential drug candidates. Due to the difficulty in obtaining skin tissue from individuals infected with EV-A71, we used hiPSC-SOs as a new strategy for studying EV-A71 virus replication and intervention targets. Combining proteomics analysis and scRNA-seq data, we revealed that EV-A71 infected multiple cell types in human skin, including epidermal cells, nerve cells, fibroblasts, and hair follicle cells and the time-dependent molecular events of EV-A71 infection were outlined. As a major cause of skin aging, senescent fibroblasts can promote skin aging by inducing chronic inflammation, reducing proliferation, and enhancing ECM degradation through enhanced activation of protein hydrolytic enzymes, including matrix-degrading metalloproteases^[180]30. By analyzing the different cell types in the scRNA-seq data, not only did we find a decrease in the number of reticular fibroblasts, we also found that collagen metabolism and inflammatory pathways (MMP1, MMP3, and MMP13) were activated in reticular fibroblasts. Meanwhile, significant downregulation of proteins associated with collagen fiber (COL1A1, COL1A2, COL2A1, COL3A1, COL5A1, COL11A1, COL14A1, and LUM) and elastin fiber-associated extracellular matrix proteins (LOX and EMILIN1) was also demonstrated in the proteomics data, which further support the results that fibroblast dysfunction leads to the decreased extracellular collagen. We concluded that EV-A71 infection may be associated with skin aging by affecting fibroblasts, which is also supported by previous studies of fibroblast senescence induced by parvovirus B19^[181]39 and influenza virus^[182]40. In addition, in the scRNA data of EV-A71-infected organoids, we also found the numbers of basal stem cells and mature EpCs decreased. Meanwhile, we observed that proteins involved in epidermal development, such as basal stem cell markers (KRT5 and KRT15) and epithelium markers (IVL and KRT12) were significantly dysregulated. These results collectively indicated a decline in the function of epidermal stem cells after viral infection, further supporting the association between viral infection and skin aging. In a previous study on EV-A71-infected human skin and oral mucosa organotypic cultures derived from prepuce and lip biopsies, viral antigens/RNA were found in the cytoplasm of epidermal and mucosal squamous cells as early as 2 dpi. In EV-A71-infected hamster models, researchers found multiple focal inflammatory skin lesions on the paws, lips, oral cavity and other parts of the skin following oral infection. Moreover, viral antigens/RNA were detected in the cytoplasm of epidermal squamous cells in these areas as well as epidermal squamous cells in oral and esophageal mucosa^[183]41. Similarly, the SCARB2 transgenic mouse model also demonstrated EV-A71 antigens in oral squamous epithelium and distal limb epidermis^[184]42. In our study, we found that EV-A71 can infect the epidermis leading to a decrease in mature epidermal cells and basal stem cells as well as a paradoxical increase in PCs, suggesting the presence of atypical cell proliferation in the epidermis. This finding suggested that EV-A71 infection may be associated with tumorigenesis. This portion of the abnormally increased PCs we identified highly expressed NNMT and SFRP1, both of which are strongly associated with cancer progression^[185]31,[186]43. The significant increase in the number of PCs was confirmed by the subsequent finding of intracellular activation of YAP/TAZ signaling which is prevalent in human malignancies^[187]32. Integrins have been reported to be important receptors for the invasion of cells by a variety of viruses such as zika virus^[188]44, Foot-and-mouth disease virus^[189]45, and SARS-CoV-2^[190]46 and this mechanotransduction is probably the main way EV-A71 activates Hippo-YAP/TAZ signaling to promote PC proliferation. As a leading cause of global cancer, certain types of viral infections can disrupt conserved signaling pathways that control cell cycle progression and promote uncontrolled cell proliferation leading to tumorigenesis^[191]47. For example, in skin tissue, human papillomavirus and herpesvirus 8 can lead to squamous cell carcinoma and Kaposi’s sarcoma, respectively^[192]48. A previous study reported that cancerous tissues had more positive EV-A71 antigen expression than adjacent non-cancerous tissues, also revealing a correlation between the virus and tumors^[193]49. Interestingly, another study reported that NNMT was overexpressed in basal cell carcinoma (BCC) of the skin and downregulated in squamous cell carcinoma (SCC). Whether this suggested that EV-A71 infection of the skin is more likely to be associated with the development of BCC rather than SCC requires further study. In addition, we also found that hair follicle stem cells are the target cells for EV-A71 infection, which may suggest that EV-A71 affects hair growth, consistent with previous reports^[194]50. Previous studies mainly used 2D cell models to investigate virus‒host interactions and antiviral targets. Here, we first utilized an organoid model infected with EV-A71 virus identification of potential drug targets. HMGB1, which is involved in both autophagy and inflammation, was identified as a universal antiviral target in EV-A71 skin infection. Inhibition of HMGB1 using a small molecular inhibitor or genetic knockdown of HMGB1 both provided favorable antiviral effects through the inhibition of virus-induced autophagy in vivo and in vitro. The inhibitor of HMGB1 (glycyrrhizic acid) is approved in Japan and China as a hepatoprotective drug for chronic hepatitis^[195]51. It has been used in a vast variety of formulations, and was reported to be anti-inflammatory, anti-ulcer, anti-allergic, antioxidant, anti-tumor, anti-diabetic and hepatoprotective. Due to these properties, its indications have been premenstrual syndrome, viral infections, hyperlipids and hyperglycemia^[196]52. The reported side effects of NSC 167409 are generally mild, including hypertension, hypokalemic and anorexia nervosa^[197]53. There was no data on its biodistribution and availability in the skin currently. But it was reported that the maximum concentrations of plasma (Cmax) of NSC 167409 after intravenous administration of 80 mg was 29.25 μg/ml^[198]54. Here, we found that the IC50 of NSC 167409 is 1.02 ± 0.10 μM, which is 34.84 times lower that its concentration in the plasma. Considering that the plasma distribution of NSC 167409 overwhelms the effective antiviral concentration, it is believed that NSC 167409 may play a role in future clinical management of enterovirus infections. It is worth noticing that NSC 167409 is also approved for children >4 years old in China^[199]55, who is permissive for enterovirus infections. Moreover, the drug reversed the abnormal proliferation of PCs in virus-infected organoids and HMGB1-mediated inflammatory responses were also inhibited. These results suggested that inhibiting HMGB1 not only hinders EV-A71 replication but also reduces atypical epidermal proliferation and skin inflammation driven by the virus. It has been reported that HMGB1 can promote tissue tumor progression through activation of YAP^[200]36,[201]56. Here, we found that inhibition of HMGB1 resulted in downregulation of YAP signaling, which also suggested a potential promotional role of HMGB1 in tumorigenesis in addition to its involvement in autophagy and inflammation (Fig. [202]6D). These results emphasize the importance of using skin organoid models as tool to aid the development of antiviral drugs, as they can provide information on crucial primary and/or secondary viral replication sites that may play an important role in viremia and person-to-person transmission. A limitation of the study is that the skin organoids in this study lack blood vessels and immune cells such as Langerhans cells and T cells which play an important role in EV-A71 infection and pathogenesis^[203]57. The goal for the future is to develop organoids with functional vascular and immune systems that will be more useful in helping us study skin infectious diseases. Methods Cells, viruses, and main chemical compounds and reagents This study was approved in accordance with the ethical standards of the institutional board of the Peking Union Medical College Hospital. Human embryonal RD cells and Vero cells were purchased from the American Type Culture Collection (ATCC) and cultured in Dulbecco’s modified Eagle’s medium (DMEM, Gibco) containing 10% fetal bovine serum (FBS, Gibco), 100 U/mL penicillin, and 100 μg/mL streptomycin at 37 °C in the presence of 5% CO[2]. EV-A71 H strains were purchased from the American Type Culture Collection (ATCC, VR-1432), and they were propagated and titrated in RD cells in DMEM containing 10% FBS^[204]58. To prepare the EV-A71 virus, the infected cell culture medium was thawed three times to release the virus and centrifuged at 3000 × g for 20 minutes to remove cellular debris. After that, the virus was purified by centrifugation at 80,000 × g for 4 hours, and then resuspended in DMEM. In the infection experiment, skin organoids were infected with 2.5 MOI of EV-A71. In the antiviral assay to evaluate the inhibition efficacy of NSC167409, the skin organoids were infected with an MOI of 10. Glycyrrhizin (NSC 167409), ammonium glycyrrhizinate, P7C3, GMX1778 (CHS828), GNE-617, OT82 and STF-118804 were purchased from Selleck Chemicals. 7-Deaza-2′-C-acetylene-adenosine (NITD008) was acquired from MedChemExpress (Shanghai, China). These reagents were dissolved in DMSO to make 100 mM stock solutions and stored at -20 °C for further study. Monoclonal mouse anti-EV-A71 virus VP1 antibody (12D7) was kindly provided by Tong Cheng Professor at Xiamen University^[205]59. Cytopathic effect inhibition assay To determine the antiviral activity of NSC 167409, a cytopathic effect (CPE) inhibition assay was performed using Vero and RD cells in a 96-well plate as previously described^[206]60. Briefly, RD cells were seeded in 96-well plates and incubated for 24 h and then treated with threefold serially diluted compounds and the virus (100 TCID[50]/ml). After incubation at 37 °C for 72 h, the cell viability of each well was measured using a CellTiter-Glo cell viability assay kit. For the cytotoxicity assay, the experimental procedure was similar to that of the CPE inhibition assay, except that the virus inoculum was replaced with DMEM containing 2% FBS. The half-maximal effective concentration (EC[50]) and half-cytotoxic concentration (CC[50]) were calculated using Origin 9.0 software. HMGB1 knockdown and western blot analysis To generate HMGB1 knockdown cells, RD cells were transfected with HMGB1-siRNA (siRNA-1: GCAGAUGACAAGCAGCCUUTT, siRNA-2: CAGGAGGAAUACUGAACAUTT) using Lipofectamine 3000 according to the manufacturer’s instructions. After transfection for 48 h, the cells were incubated with the EV-A71 virus at an MOI of 0.1 PFU/mL and then harvested after 24 h of infection. To verify the knockout efficiency, a western blot analysis was carried out. The total proteins were separated on a 12% SDS‒PAGE gel and then transferred onto a PVDF membrane. After being blocked for 1 h, the membranes were incubated with anti-HMGB1, EV-A71-VP1, and LC3B antibodies or an anti-α-tubulin antibody as a control at 4 °C overnight. The membranes were washed three times and incubated with the appropriate secondary antibody for 1 h at room temperature. Finally, the protein bands were imaged using a FluorChem imaging system. Generation of skin organoids The hiPSC line (passages 12-35) was purchased from Shownin Biotechnologies Co., Ltd (RC01001-B, Female). Cells were cultured on 6-well plates coated with Corning® Matrigel® using mTeSR™ Plus medium. The medium was replenished every other day with 3 ml each time. When reaching 70% ~ 80% confluency, the cells were passaged in tiny clusters with gentle cell dissociation buffer. Dissociated cell colonies were passaged in mTeSR™ Plus medium containing 10 μM Y27632. Twenty-four hours after passage, medium was replenished with fresh mTeSR™ Plus medium without Y27632. Skin organoids were generated as described in previous research^[207]18. Briefly, the hiPSCs were dissociated into single cells and seeded into low attachment 96-well U-bottom plates at 3500 cells per well with 100 μl Essential 8 (E8) + 20 μM Y27632 to form aggregates. After 24 h, the aggregates were supplemented with 100 μl E8 for another day of incubation. On day 0 of differentiation, the aggregates were transferred to a new 96-well plate in 100 μl Essential 6 (E6) containing 2% Matrigel, 10 μM SB-431542, 4 ng/mL bFGF and 2.5 ng/mL BMP4 for 3 days. Subsequently, the aggregates were treated with 25 µl of E6 containing 250 ng/mL bFGF and 1 μM LDN for another 3 days. On day 6, 75 μl of fresh E6 was added to provide nutrition, and half of the medium was replaced with fresh E6 every 2 days until day 12. On day 12, the aggregates were transferred to a low attachment 24-well plate in 500 μl OMM (50 ml) composed of 24.5 ml advanced DMEM-F12, 24.5 ml neurobasal medium, 500 μl B-27 supplement without vitamin A, 500 μl Glutamax supplement, 250 μl N2 supplement, 91 μl 2-mercaptoethanol, and 100 μl Normocin with 1% Matrigel and placed on an orbital shaker at 65 rpm. On day 15, half of the medium was replaced with fresh OMM with 1% Matrigel. From day 18 to organoid maturation, half of the medium was replaced with fresh OMM every 2-3 days. On day 75, TGFBI (Transforming Growth Factor-Beta-Induced Protein Ig-H3) was added in the culture medium to promote epidermal development. The mature skin organoids cultured for 143 days were used for all experiments in this study, including the EV-A71 infection, scRNA-seq, and proteomics analysis. hiPSC-derived skin xenograft infection All animal work was approved by the Institutional Animal Care and Use Committee of the Beijing Institute of Pharmacology and Toxicology (approval number: IACUC-2022-039W). A total of 15 male 6-week-old BALB/c nude mice were used in this study. The mice were anesthetized with sodium pentobarbital (1% v/v diluted in normal saline). Afterwards, a small incision was made on the back skin of each mouse and two skin organoids were subcutaneously transplanted into each body. The mice were maintained at 22 °C with 40–60% humidity and a 12 light/12 dark cycle and the postoperative status of them was observed daily. One month later, the mice were further used for drug evaluation. Briefly, 3 groups were set, termed control group, EV-A71 infection group, and inhibitor-treatment group (n = 5 for each group). In treatment group, the prophylactic treatment was performed on mice with NSC167409 (2 mg/kg) by subcutaneous injection. Four hours after drug administration, the EV-A71 virus was directly inoculated into the xenograft at 7×10^6 PFU per mouse except for the control group. Histology and immunofluorescence The samples were fixed in 4% formaldehyde at 4 °C overnight and dehydrated using an alcohol concentration gradient. After xylene treatment and paraffin-embedding, samples were cut into 4-μm-thick sections for immunofluorescence analysis. For general immunofluorescence, the sections were deparaffinized in xylene and rehydrated using an alcohol gradient. Antigen retrieval was performed by boiling the sections with citrate buffer in a microwave for at least 12 min. Goat serum (10%) was added and incubated for 1 h at room temperature, followed by overnight incubation with primary antibodies at 4 °C. The next day, the sections were incubated with secondary antibodies for 1 h at room temperature and then stained with DAPI. After sealing with antifade mounting medium, images were acquired at different magnifications and analyzed using a Nikon A1R + N-STORM confocal microscope. Multiplex immunofluorescence was performed using tyramide signal amplification (TSA)-dendron-fluorophores and a NEON Allround Discovery Kit for FFPE (Histova Biotechnology). Briefly, the sections were treated with 3% H[2]O[2] for 10 min at room temperature after antigen retrieval to remove endogenous peroxidases. The primary antibody was added and incubated overnight at 4 °C or 2 to 3 h at 37 °C, followed by incubation with HRP-conjugated secondary antibody and subsequent TSA-dendron-fluorophore signal amplification. The sections were boiled in elution buffer for 15 s of to degrade the primary antibody and start a new round of antibody incubation. Different antigen targets were stained with different TSA fluorescence through multiple rounds of TSA signal amplification. Multiplex immunofluorescence images were acquired using a Zeiss LSM880 + ELYRAS.1 confocal microscope. Proteomics sample preparation The EV-A71-infected hiPSC-SOs at 2, 4, and 8 dpi and hiPSC-SOs without infection were collected to for proteomics analysis. Samples from the skin organoids were scraped into a new EP tube with 20 μL of urea buffer (8 M urea, 150 mM Tris-HCl, 10 mM DTT, pH 8.0). An additional 10 μL of buffer was added to the tube. Steel balls were added to 30 μL of buffer to homogenize (70 Hz) the tissue for 1 min. After centrifugation at 14,000 × g for 10 min at 4 °C, the supernatants were transferred to clean tubes. Next, the extracted proteins were reduced at 37 °C for 1 h and alkylated in 25 mM iodoacetamide at room temperature for 30 min in the dark. The protein samples were digested with Lys-C (1 μg at 37 °C for 4 h) and trypsin (enzyme-to-substrate ratio of 1:50 at 37 °C for 16 h), desalted using C18 cartridges, and vacuum-dried using a SpeedVac. The peptides from each sample were collected for subsequently MS analyses, while a part of the peptides from all the samples were pooled and fractionated for transition library construction. The pooled peptides mixtures were fractionated using a C18 column (Waters BEH C18, 4.6×250 mm, 5 µm) on a Rigol L3000 HPLC instrument operating at 1 mL/min, with the column oven set to 50°C. Mobile phases A [2% acetonitrile (ACN), adjusted pH 10.0 using ammonium hydroxide] and B (98% ACN, adjusted pH 10.0 using ammonium hydroxide) were used to develop a gradient elution. The solvent gradient was set as follows: 5% B, 0 min; 5–8% B, 5 min; 8–18% B, 35 min; 18–32% B, 22 min; 32–95% B, 2 min; 90% B, 4 min; 95–5% B, 4 min. The eluates were monitored at UV 214 nm, collected in a tube, and merged into three fractions. The peptide fractions were dried in a vacuum. Finally, the peptides of each sample and the three peptide fractions of pooled samples were reconstituted in 0.1% (v/v) formic acid in water. Then, 0.2 µL of standard peptides (iRTkit, Biognosys) was added to the peptide sample for subsequent analyses. Mass spectrometry analysis The peptides from all the samples were pooled for transition library construction with shotgun proteomics analyses were performed using an EASY-nLCTM 1200 UHPLC system coupled with an Orbitrap Q Exactive HF-X mass spectrometer (Thermo Fisher Scientific) operating in the data-dependent acquisition (DDA) mode. The data-independent acquisition (DIA) scan mode was used for single-shot for all samples. The fractionated samples of the pool were acquired using the top 40 DDA scan mode. Both acquisition schemes were combined with the same liquid chromatography gradient. The mass spectrometer was operated using Xcalibur software (version 4.1). The DDA scan settings on the full MS level included an ion target value of 3 × 10^6 charges in the 350-1,500 m/z range, with a maximum injection time of 80 ms and a resolution of 120,000 at 200 m/z. At the MS/MS level, the target value was 5 ×10^4 charges with a maximum injection time of 45 ms and a resolution of 15,000 at m/z 200. For the MS/MS events only, precursor ions with 2-7 charges not on the 16 s dynamic exclusion lists were isolated within a 1.6 m/z window. Fragmentation was performed by higher-energy C-trap dissociation with normalized collision energy of 27 eV. DIA was performed with one full MS event, followed by 42 MS/MS windows in one cycle. The full MS settings included an ion target value of 3 × 10^6 charges in the 350-1,500 m/z range, with a maximum injection time of 50 ms and a resolution of 60,000 at m/z 200. The DIA precursor windows ranged from 378 m/z (lower boundary of the first window) to 1345 m/z (upper boundary of the 42nd window). The MS/MS settings included an ion target value of 1 × 10^6 charges for the precursor window, with an Xcalibur-automated maximum injection time and a resolution of 30,000 at m/z 200. Proteomic MS/MS data processing The MS data of the fractionated pool (DDA data, three fractions) were analyzed by Proteome Discoverer (version 2.4). Then, the DDA data and the single-shot samples (DIA data) were used to generate a DDA library and direct DIA library, respectively, which were computationally merged into a hybrid library in Spectronaut (version 14.9.201124.47784, Biognosys). Then, the raw DIA data were processed on Spectronaut using the default settings. The MS/MS data files generated from human skin organoids were searched against the UniProtKB reviewed human protein sequences (Swiss-Prot) (downloaded in Sep. 2022, containing 20,398 entries) combined with the reviewed and unreviewed EV-A71 protein sequences in UniProtKB (downloaded in Jul. 2023, containing 14,209 entries). The searches used carbamidomethylation as the fixed modification and acetylation of the protein N-terminus and oxidation of methionines as the variable modifications. Default settings were used for other parameters. In brief, a trypsin/P proteolytic cleavage rule was used, permitting a maximum of two miscleavages and a peptide length of 7–52 amino acids. Protein intensities were normalized using the “Local Normalization” algorithm in Spectronaut based on a local regression model; the retention time prediction type was set to dynamic iRT and correction factor for a window. The mass calibration was set to local mass calibration. Decoy generation was set to Inverse. Interference correction on the MS2 level was enabled, removing fragments for quantification based on the presence of interfering signals but maintaining at least three fragments for quantification. The false discovery rate was estimated using the mProphet approach and set to 1% at the peptide level. Skin organoid cell collection, scRNA-seq, and scRNA-seq data analysis The EV-A71-infected hiPSC-SOs at 8 dpi and hiPSC-SOs without infection were collected to generate single-cell gene expression libraries. The collected organoids were first washed with ice-cold RPMI 1640 and then incubated with TrypLE™ express enzyme at 37 °C with gentle swirling until the organoids dissociated into a single-cell suspension with no obvious cell aggregation. After neutralizing the trypsin with bovine serum albumin (BSA) solution, the suspension was filtered through a 40 μm Flowmi cell strainer to remove cell debris. Then, the cell count and viability were estimated using a fluorescence Cell Analyzer (Countstar® Rigel S2) with acridine orange/propidium iodide reagent to ensure that the cell concentration and viability were suitable for sequencing. The fresh cells were washed twice with RPMI 1640 and resuspended at 1×10^6 cells per ml in 1×PBS and 0.04% BSA. Finally, the single-cell RNA-Seq libraries were prepared using a SeekOne® MM Single Cell 3’ library preparation kit (SeekGene Catalog No. [208]K00104) and then sequenced on an Illumina NovaSeq 6000 with PE150 read length. Cell Ranger software (10× Genomics) ([209]https://support.10xgenomics.com/v6.1.1) was used to process the raw scRNA-seq data. The FASTQ files were aligned to the human genome and transcriptome (hg38) to generate a gene expression matrix. The R package Seurat (v.4.3.0) was used for further data analysis. First, cells with less than 500 or more than 5000 detected genes were filtered, and those with more than 15% mitochondrial detected genes were further excluded. After quality control, a total of 23,483 cells remained. After deleting unnecessary cells, the “NormalizedData” function was used to log-normalize the count assay, and the top 2000 highly variable genes were determined using the “FindVariableFeatures” function for downstream bioinformatics analyses. Next, we continued processing the dataset by applying the linear transformation “ScaleData” function. The batch effects among samples were removed using the “RunHarmony” function. The “FindAllMarkers” function was used to identify cluster markers with the following settings: min percentage of the gene expressed = 0.1, log fold change threshold = 0.25 and Bonferroni-adjusted p < 0.05. Differentially expressed genes (DEGs) between the two groups of cells were identified from the scRNA-seq data using the two-sided Wilcoxon rank-sum test. “FindMarkers” function was used to find the upregulated or downregulated genes in epidermal cells and fibroblasts in control and EV-A71-infected hiPSC-SOs. The genes with Bonferroni-adjusted p value less than 0.05 and |logFold Change | > 0.25 were identified as DEGs. Statistical analysis The protein quantification values were normalized by multiplying the fraction of the total intensity by 10^6 followed by log2 transformation. Pairwise comparisons to determine the differentially expressed proteins (DEPs) between the hiPSC-SOs with EV-A71 infection at 2 (n = 3), 4 (n = 3), and 8 (n = 3) dpi and the hiPSC-SOs without EV-A71 infection (control group, n = 3) were performed by the two-sided moderated t test using the R package Limma (version 3.50.3). The significantly DEPs were determined according to the Benjamin–Hochberg (BH) adjusted p value of <0.01 and were considered statistically significant, and those with log2 EV-A71 infection/Control > 0.585 or < -0.585 were considered upregulated or downregulated, respectively. One-way ANOVA was used to determine whether there were statistically significant differences in the normalized intensities among proteins identified in hiPSC-SOs with different EV-A71 infection times. The specifically highly and lowly expressed proteins at different infection timepoint were defined as those with BH adjusted one-way ANOVA p < 0.01 and log2FC > 0.585 or < -0.585, respectively. The fold change (FC) represents the ratio of the normalized intensity of the protein identified in the hiPSC-SOs at the indicated infection time to those at all other samples. For other experimental results, the statistical data are shown as the mean ± standard deviation (SD). The unpaired two-tailed Student’s t test was performed to compare data between groups. Each experiment was performed at least three times (n ≥ 3). The results were considered significant at p ≤ 0.05 and the exact p values are provided in Source Data. Statistical tests were conducted using GraphPad Prism (version 9.0.0). Principal coordinates analysis (PCoA) of proteins was performed using the R package ape (version 5.5), and the values of each sample were validated^[210]61. The volcano plot and heatmap of significant differentially expressed proteins were constructed using the R package ggplot2^[211]62 (version 3.3.5) and ComplexHeatmap^[212]63 (version 2.8.0, distance: Pearson, linkage: complete). The online tool [213]DAVID^[214]64 was used for function enrichment analysis of the DEGs and DEPs via biological process of Gene Ontology^[215]65, and the Kyoto Encyclopedia of Genes and Genomes^[216]66 pathway enrichment analysis was used to identify the pathways that the DEPs were associated with ([217]https://www.kegg.jp/kegg/pathway.html). The statistical significance of enrichment analysis is performed using the one-sided hypergeometric test and the exact p values are provided in Source Data. Reporting summary Further information on research design is available in the [218]Nature Portfolio Reporting Summary linked to this article. Supplementary information [219]Supplementary Information^ (8.6MB, pdf) [220]41467_2025_57610_MOESM2_ESM.docx^ (13.4KB, docx) Description of Additional Supplementary Files [221]Supplementary Data 1^ (535KB, xlsx) [222]Supplementary Data 2^ (282.1KB, xlsx) [223]Supplementary Data 3^ (211.5KB, xlsx) [224]Supplementary Data 4^ (350KB, xlsx) [225]Supplementary Data 5^ (19KB, xlsx) [226]Reporting summary^ (8.2MB, pdf) [227]Transparent Peer Review file^ (1.2MB, pdf) Source data [228]Source data^ (25.2MB, xlsx) Acknowledgements