Human skeletal muscle organoids model fetal myogenesis and sustain uncommitted PAX7 myogenic progenitors

In vitro culture systems that structurally model human myogenesis and promote PAX7+ myogenic progenitor maturation have not been established. Here we report that human skeletal muscle organoids can be differentiated from induced pluripotent stem cell lines to contain paraxial mesoderm and neuromesodermal progenitors and develop into organized structures reassembling neural plate border and dermomyotome. Culture conditions instigate neural lineage arrest and promote fetal hypaxial myogenesis toward limb axial anatomical identity, with generation of sustainable uncommitted PAX7 myogenic progenitors and fibroadipogenic (PDGFRa+) progenitor populations equivalent to those from the second trimester of human gestation. Single-cell comparison to human fetal and adult myogenic progenitor /satellite cells reveals distinct molecular signatures for non-dividing myogenic progenitors in activated (CD44High/CD98+/MYOD1+) and dormant (PAX7High/FBN1High/SPRY1High) states. Our approach provides a robust 3D in vitro developmental system for investigating muscle tissue morphogenesis and homeostasis.


Introduction
Novel skeletal muscle model systems are required to further elucidate the process of human myogenesis as well as investigate muscular disorders and potential gene, cell, or drug therapies.Twodimensional (2D) culture conditions guide pluripotent stem cell (PSC) differentiation toward skeletal muscle lineage using sequential growth factor applications and/or conditional PAX7 expression (Chal et al., 2015;Xi et al., 2017;Shelton et al., 2014;Borchin et al., 2013;Darabi et al., 2012).Further, surface marker expression can be utilized to isolate myogenic progenitors with in vivo repopulation potential (Magli et al., 2017;Hicks et al., 2018;Al Tanoury et al., 2020;Sun et al., 2022).While the few described three-dimensional (3D) differentiation approaches have provided cohorts of terminally differentiated myofibers, focus on potential interactions with the vasculature and nervous system has neglected assessment of the developmental identity or sustainability of myogenic progenitors (Faustino Martins et al., 2020;Maffioletti et al., 2018;Rao et al., 2018).Single-cell technologies increasingly provide databases for deciphering myogenic trajectories and expression profiles of myogenic stem and progenitor cells (Barruet et al., 2020;Rubenstein et al., 2020;Xi et al., 2020), enabling full evaluation of the ability of PSC differentiation protocols to mimic human development.Translation to model muscular dystrophies and investigate potential interventions in vitro necessitates methods that provide expandable populations of muscle progenitors while promoting self-renewal and preserving a quiescent, non-dividing, state (Quarta et al., 2016;Montarras et al., 2005).
Most vertebrate skeletal muscle progenitors develop from the paraxial mesoderm via transient embryonic developmental structures (somites and dermomyotome) into the skeletal muscle system that spans the whole body.Here, we evaluate human skeletal muscle organoids as a novel system to structurally investigate myogenic differentiation from human induced pluripotent stem cells (iPSCs) in a 3D environment, mimicking pathways described for chicken and mouse (Buckingham and Rigby, 2014).We develop a comprehensive supplementation/reduction protocol to first drive differentiation toward paraxial mesoderm through application of the GSK3 inhibitor CHIR99021, BMP inhibitor LDN193189, and bFGF.Subsequent stimulation with WNT1A, SHH, FGF, and HGF is designed to promote derivation of organized structures reassembling neural plate border and dermomyotome.We then aim to arrest neural lineage via FGF removal, while stimulating with HGF to selectively promote propagation of undifferentiated myogenic progenitors and consequent generation of fetal myofibers.Our goal is to provide PAX7 + myogenic progenitors in a non-dividing quiescent state sustainable over weeks after differentiation induction.Single-cell analysis will position cells along the quiescentactivation myogenic trajectory discriminating dormant (PAX7+, FBN1+, SPRY11+, CHODL1+) and activated (CD44+, CD98+, MYOD1+, VEGFA+) states.We thus seek to develop and validate a new skeletal muscle organoid system for investigating human myogenesis with translational potential for disease modeling and therapy development.
Consequently, for favoring myogenesis we stimulated organoid culture with FGF and HGF (Chargé and Rudnicki, 2004;Figure 1A).Surprisingly, on day 17, organoids constituted a mosaic of neural crest and myogenic progenitor cells.Cells with epithelial morphology were TFAP2A + , SOX2 + , PAX3 + , and PAX7 + , indicating formation of neural plate border epithelium (Roellig et al., 2017;Figure 1G).In cell clusters with mesenchymal morphology, we detected specified migrating PAX3 + /SOX10 + neural crest progenitors and PAX3 + /SOX2 -/SOX10 -cells of myogenic origin (Figure 1G).At this stage, eLife digest Humans contains around 650 skeletal muscles which allow the body to move around and maintain its posture.Skeletal muscles are made up of individual cells that bundle together into highly organized structures.If this group of muscles fail to develop correctly in the embryo and/or fetus, this can lead to muscular disorders that can make it painful and difficult to move.
One way to better understand how skeletal muscles are formed, and how this process can go wrong, is to grow them in the laboratory.This can be achieved using induced pluripotent stem cells (iPSCs), human adult cells that have been 'reprogrammed' to behave like cells in the embryo that can develop in to almost any cell in the body.The iPSCs can then be converted into specific cell types in the laboratory, including the cells that make up skeletal muscle.
Here, Mavrommatis et al. created a protocol for developing iPSCs into three-dimensional organoids which resemble how cells of the skeletal muscle look and arrange themselves in the fetus.To form the skeletal muscle organoid, Mavrommatis et al. treated iPSCs that were growing in a threedimensional environment with various factors that are found early on in development.This caused the iPSCs to organize themselves in to embryonic and fetal structures that will eventually give rise to the parts of the body that contain skeletal muscle, such as the limbs.Within the organoid were cells that produced Pax7, a protein commonly found in myogenic progenitors that specifically mature into skeletal muscle cells in the fetus.Pax 7 is also present in 'satellite cells' that help to regrow damaged skeletal muscle in adults.Indeed, Mavrommatis et al. found that the myogenic progenitors produced by the organoid were able to regenerate muscle when transplanted in to adult mice.
These findings suggest that this organoid protocol can generate cells that will give rise to skeletal muscle.In the future, these lab-grown progenitors could potentially be created from cells isolated from patients and used to repair muscle injuries.The organoid model could also provide new insights in to how skeletal muscles develop in the fetus, and how genetic mutations linked with muscular disorders disrupt this process.myogenic lineages appeared to be primarily represented by PAX3 + (9.35 ± 0.07%) rather than PAX7 + cells, as the PAX7 expression pattern (15.3 ± 0.1% PAX3 + /PAX7 + , 5.11 ± 0.13% PAX7 + ) predominantly overlapped with that of SOX2 and TFAP2A (Figure 1G and Figure 1-figure supplement 1B).Morphologically, neural plate border and dermomyotomal populations exhibited uneven distribution, and thereby subsequent neural crest and myogenic migration processes introduced organoid size variability (Figure 1B and G).From day 17 onward, we omitted FGF stimulation to cease neural crest development and promote delamination/migration of PAX3 + /PAX7 + /SOX10 -progenitor cells (Murphy et al., 1994;Figure 2A).Strikingly, until day 23, we observed committed myogenic populations through detection of fast MyHC myofibers in the proximity of PAX7 + but SOX2 -/TFAP2A -/SOX10 - negative cells (Figure 1H, Figure 1-figure supplement 1C-E).Consistently, expression profiling indicated downregulation of neural tube/crest lineage markers and significant upregulation of muscle precursor migrating markers (Buckingham and Rigby, 2014), such as LBX1, CXCR4, and MEOX2 (Figures 1E and 2A,. At 8 wk, organoids showed profound changes in transcription profiling (Figure 2, Figure 2-figure supplement 1B).Gene ontology enrichment analysis highlighted an ongoing developmental transition with muscle development, muscle contraction, and myofibril assembly among the top upregulated gene ontology terms and neurogenesis and nervous system development among the top downregulated terms (Figure 2B).In addition, we detected downregulation of key markers characterizing neural tube and neural crest development (Soldatov et al., 2019), such as PAK3, DLX5, B3GAT1, FOXB1, HES5, and OLIG3 (Figure 2, Figure 2-figure supplement 1C).Interestingly, we could spatially visualize this process using immunocytochemistry for SOX2, TFAP2A, and SOX10-expressing cells that were restricted to the inner portion of the organoid, and probably not susceptible to culture stimulation at 5 wk (Figures 1H,2C and D ,.This neural/myogenic lineage spatial orientation could be visualized even on day 84 through the presence of TUJ1 + neurons confined to inner organoid areas and close to SOX2 + epithelium, while FastMyHC + myofibers occupied exterior areas (Figure 2E).On the other hand, substantial proportions of migratory cells in proximity of FastMyHC + , MF20 + myofibers expressed PAX7 (40 ± 0.3%) but not SOX2, TFAP2A, SOX10, or MYOD1 (3.55 ± 0.32%) (Figure 2F, Figure 2-figure supplement 1D and E).This behavior is further illustrated by the presence of MYOD1 + cells confined to organoid periphery (Figure 2D).

Functionality and maturation of organoid-derived myofibers
Regarding localization and functionality of organoid-derived myofibers, immunocytochemistry revealed positive staining for dystrophin and a continuous laminin sheath around individual muscle fibers.Ultrastructure and two-photon microscopy analysis depicted well-developed sarcomeres (Figure 3, Figure 3-figure supplement 3A-D).Moreover, patch-clamp experiments, upon superfusion of cells with acetylcholine (ACh), indicated inward current generations that rapidly declined to a low steady-state level (Figure 3K, Figure 3-figure supplement 3E).The I/V curve of this current component showed almost linear behavior and reversal potential around 0 mV (Figure 3L).These data are in line with previous studies that described currents with analogous properties as characteristic nAChR currents (Jahn et al., 2001;Shao et al., 1998).Application of a fluorescent biosensor for monitoring changes in cytosolic Ca 2+ (Twitch2B) revealed nAChR efficiency in modulating intracellular expression of FastMyHC + and PAX7 + myogenic populations, while SOX2 + neural populations demarcate SOX2 neural plate border epithelium location as observed at earlier stages (day 16) (C); PAX7 cells are of myogenic origin (PAX7 + /SOX10 -), MF20 + myotubes are in their proximity and MYOD1 + cells appear at organoid periphery (D); TUJ1 + neurons are restricted to inner organoid areas and close to SOX2 + epithelium, while FastMyHC + myofibers occupy exterior organoid areas (E).(F) Histographs based on FACS intracellular quantification depicting percentage of PAX7 + or MYOD1 + cells through differentiation protocol.For each replicate, 10 organoids were pooled, n = 10.Statistics: *p<0.05,**p<0.01,***p<0.001,****p<0.0001,ns: not significant.(G) Heatmap of HOX gene cluster emphasizes organoid culture limb axial anatomical identity by depicting transition from an initial limb bud (HOX 9-10) toward a more distal identity   at 8 and 16    ), following application of ACh (10 μM), illustrated a rapid increase in FRET ratio that gradually declined in the presence of ACh, probably reflecting desensitization of nAChRs.These results demonstrated that nAChR in skeletal muscle cells are functional in terms of inducing Ca 2+ release from intracellular stores.
Using our organoid protocol, we successfully derived fetal muscle progenitors and electrophysiologically functional myofibers from hiPSC lines with wild type and Duchenne muscular dystrophy genetic backgrounds (Figure 3, Figure 3-figure supplement 3E and F).

Organoid-derived myogenic progenitor comparison to human fetal and adult progenitor/stem cells
To evaluate developmental identity of the myogenic cluster, we isolated ITGβ1 + /CXCR4 + organoidderived myogenic progenitors via FACS (Garcia et al., 2018) and compared to human fetal muscle progenitors, demonstrating high similarity (Pearson'ś correlation co-efficiency, rho = 0.9), with myogenic markers such as PAX7, MYF5, MYOG, and MYOD1 at comparable expression levels (Figure 4-figure supplement 3A-F).Differential expression comparison verified expression of extracellular matrix collagens and proteins, such as COL4A1, COL5A1, COL6A1, COL15A1, FBN1, and CHODL, in myogenic progenitors similar to 17th week human fetal tissue (Pearson'ś correlation co-efficiency, rho = 0.97) (Figure 4E).Further, to evaluate myogenic potency in vitro, isolated ITGβ1 + /CXCR4 + organoid-derived myogenic progenitor cells were re-plated and allowed to differentiate under the same organoid culture protocol, which demonstrated capacity to generate Fast-MyHC+ and PAX7+ populations within 10 d (Figure 4-figure supplement 3B and C).Subsequently, absence from FastMyHC + positive myofibers is highlighted.(C) FBN1 + microfibrils are located toward organoid interior.(D) Pseudotime ordering for myogenic progenitors and myoblast corresponding clusters highlights distinct developmental trajectories promoting myogenic commitment and self-renewal.(E) Correlation coefficient plot for Log2 fold change (Log2 FC) values for isolated myogenic progenitors from human fetal tissue (17w) and FACS-sorted CXCR4 + /ITGβ1 + organoid-derived myogenic progenitors (16w).PAX7, COL1A1, COL3A1, COL4A1, COL5A1, COL15A1, FBN1, and CHODL and further extracellular matrix-related genes are highlighted on the plot.Pearson's correlation coefficient, rho = 0.9 for global expression comparison and rho = 0.97 for selected genes.(F) UMAP color-based clustering divides non-dividing myogenic progenitors and adult satellite cells into four clusters with distinct molecular signatures: satellite cells characterized by SPRY1 high /ZFP36 high /MYF5 high expression, co-clustered with dormant SPRY1 high /FBN1 high /CHODL high /PAX7 high , activated CD44 high /CD98 high /MYOD1 + , and committed, NEB high /MYH3 high /MYOD1 high organoid-derived myogenic progenitors.Dots correspond to adult satellite cells from GSE130646 database, triangles correspond to organoid-derived myogenic progenitors.(G) Violin plots depicting relative expression levels of markers specific for quiescent PAX7 or activated MYF5 muscle stem cell state across adult satellite cells (GSE130646) and organoid-derived myogenic progenitors subclusters.(H) Ridge plots of developmental score distribution of myogenic progenitors across in vivo or in vitro stages based on the difference between upregulated satellite cell and embryonic markers from human reference atlases for weeks (Wk) 5-18 embryonic and fetal stages, years (Yr) 7-42 adult satellite cells and skeletal muscle (SM) organoid.(I, J) In vivo engraftment potential of human myogenic progenitors.100.000CD82+ sorted human cells were injected into Tibialis anterior muscle of nude mice.(I) Control mice were not injected (J).Six weeks post transplantation, transverse cryosections of muscle fibers were stained with huLamin A/C (green), dystrophin (red), and DAPI (blue).Human cells appear green and red in contrast to murine cells, which only show a dystrophin-positive red staining.Scale bars 200 μm in (I, J).
The online version of this article includes the following figure supplement(s) for figure 4:        comparison to available transcriptomic dataset of human adult satellite cells at single-cell resolution divided myogenic progenitors and adult satellite cells into four clusters with distinct molecular signatures (Figure 4F).Interestingly, myogenic progenitors were enriched for extracellular matrix proteins, while satellite cells mainly upregulated genes repressing transcription/translation, such as ZFP36 and TSC22D1, or related to early activation response, such as MYF5, JUN, and FOS (Figure 4-figure supplement 4A and B).In line, organoid-derived myogenic progenitors exhibited higher NOTCH signaling activity in comparison to satellite cells, with NOTCH3 and RBPJ being among enriched markers (Figure 4-figure supplement 4B).In contrast, adult satellite cells exhibited PAX7 low /MYF5 high expression profiling, presumably due to tissue dissociation, thereby indicating a tendency for activation rather than preservation or reinstatement of quiescency (Machado et al., 2017;Seale et al., 2000;Figure 4G).Pseudotime ordering showed two distinct clusters, with adult satellite cells downstream from non-dividing myogenic progenitors (Figure 4-figure supplement 4C).Consistently, downregulation of genes promoting quiescence like PAX7, NOTCH3, and RBP, and upregulation of activation genes like MYF5, JUN, and FOS along the trajectory (Figure 4figure supplement 4D), was a further indication that organoid-derived myogenic progenitors resided in dormant non-dividing state and that our organoid platform promoted sustainability of myogenic progenitors.
Cell-cell communication analysis of organoids at 12 wk indicates that myogenic progenitors influence their own fate, mainly with extracellular matrix (ECM)-related signals, as well as receive signals predominantly from the mesenchyme but not the neural progenitor cluster (Figure 4-figure supplement 5).
To evaluate for reproducibility of organoid development, we applied diffusion map analysis on qPCR-based expression analysis of 32 genes at early stages and integrative analysis on scRNAseq datasets of mature stages of organoid development.The data indicate highly conserved cluster representation of myogenic progenitors at all states, together with skeletal muscle myofibers, fibroadipogenic progenitors, and neural progenitors-related clusters (Figure 4-figure supplement 6).
In addition, we have performed a developmental score distribution of myogenic progenitors based on the difference between upregulated satellite cell and embryonic markers from the human reference myogenic atlases (Xi et al., 2020) and adult satellite cell data (Rubenstein et al., 2020;Xi et al., 2020) in comparison to our organoid protocol (Figure 4H).Organoid-derived myogenic progenitors represent the late fetal stage of maturation partially overlapping with adult satellite cell scoring.We note the heterogeneity of adult satellite cell populations when performing developmental scoring in line with recent reports (Barruet et al., 2020).
Finally, we transplanted CD82-positive progenitors from our organoids into the Tibialis anterior (TA) muscle of immunodeficient mice to complement our study with an in vivo experiment (Alexander et al., 2016;Marg et al., 2019;Al Tanoury et al., 2020).The CD82 positivity used for FACS selection of myogenic progenitors prior to transplantation almost exclusively overlaps with Pax7-positive cells, being a subcluster of them (Figure 4-figure supplement 4E).Six weeks post transplantation, we could verify clusters of human Lamin A/C-positive cells in the transplanted but not in the control group (Figure 4I and J).We have measured the size of myotubes with 41 ± 6 µm for the human and 63 ± 7 µm for the mice mean diameters (n = 15 each).

Discussion
Human skeletal muscle organoids offer a new cell culture system to study human myogenesis, in particular fetal myogenic progenitors.We demonstrate that modulation of Matrigel-embedded embryonic bodies with WNT, BMP, and FGF signaling at early stages leads to paraxial mesoderm formation (Figure 1B).Further, under guided differentiation, we could promote concomitant development of neural and paraxial mesodermal lineages and derive mesodermal populations with somitic and dermomyotomal origin (Figure 1C-F).With WNT1A and SHH stimulation, neural lineage is directed toward dorsal neural tube/crest development which benchmarks the structural recapitulation of neural plate border epithelium (Figure 1G).In vitro neural plate border can be utilized to assess transcriptomic networks and cell fate decision during human neural crest formation.
Delaminating from the dermomyotome, undifferentiated PAX3 progenitor cells reorientate and align in the cranio-caudal axis to form the underlying myotome or completely detach from the dermomyotome and migrate into the surrounding mesenchyme at the limb sites, where they propagate and become committed to skeletal muscle progenitors (Relaix et al., 2005).By stimulating organoid culture at the neural plate border/dermomyotomal stage with bFGF/HGF, we could further visualize both migration of myogenic progenitors and migration/specification of neural crest populations (Figures 1A, G, H and 2A).Further, by omitting FGF during organoid development, we could detect a continuous upregulation of genes involved in the myogenic migration process, such as LBX1, PAX3, PAX7 and CXCR4, but not in genes characterizing neural tube or neural crest development, such as SOX10, TFAP2A, PAK3, DLX5, B3GAT1, FOXB1, HES5, and OLIG3.This indicates that organoid culture conditions and specifically stimulation with HGF favored skeletal muscle over neural lineage (Figure 2C-E).Interestingly, we could show that by stimulating organoid culture with SF/HGF, an essential signal for detachment and subsequent migration, but not for specification of the cells at the lateral lip of the dermomyotome (Dietrich et al., 1999), we could preserve the PAX3+/PAX7+ myogenic migratory population in an undifferentiated and uncommitted state (Figure 2D and E).Strikingly, expression profiling based on HOX gene cluster supported this notion as over time the organoid culture adopted more distal than proximal limb axial anatomical identity (Figure 2F).
Fetal development is characterized by PAX3+ to PAX7+ myogenic progenitor transition (Seale et al., 2000;Relaix et al., 2005), which we were able to demonstrate in our organoid culture.Our data further support organoid culture representing fetal stages as we could detect NFIX upregulation and verify the presence of myofibers expressing fetal MyHC isoform as well as NCAM, M-Cad, or MCAM in the proximity of PAX7 progenitors (Figure 3G-J).Consequently, transcriptome comparison indicated high similarity to human fetal muscle tissue (17 wk, Pearson correlation, rho = 0.9) (Figure 4E), as well as expression of several satellite cell markers (Figure 4-figure supplement 1A).This was further verified by comparison to a single-cell dataset from the human skeletal muscle cell atlas (Figure 4, Figure 4-figure supplement 4) Interestingly, single-cell resolution showed adult satellite cell clusters with organoid-derived myogenic progenitors (Figure 4F).In addition, pseudotime ordering indicated that organoid-derived myogenic progenitors reside developmentally upstream of activated satellite cells, upregulating markers associated with quiescence such as NOTCH signaling and extracellular matrix proteins (Figure 4-figure supplement 2E and F).Preservation of myogenic progenitors in a non-dividing state without activating MYOD1 and self-renewal (Figure 4D, Figure 4-figure supplement 2E and F) appeared to be responsible for the observed sustainable propagation of PAX7 progenitor cells even 14 wk post differentiation (Figures 2F and 3A).
Patterning in our 3D protocol provides progenitors in a more mature late fetal state partially overlapping with adult satellite cell developmental scoring (Figure 4H).The 'uncommitted' Pax7 progenitor status is demonstrated by our precursor populations bulkRNAseq and scRNAseq profiling (Figure 4A, E and F).In this context, we could observe high expression of extracellular matrix proteins and upregulated NOTCH signaling in dormant non-dividing myogenic progenitors (Figure 4-figure supplement 2A and E).This phenotype is similarly described for human fetal tissue myogenic progenitors (Figure 4E, Pearson correlation, rho = 0.97, Figure 4-figure supplement 5).Studies evaluating engraftment potential of satellite cells and fetal muscle progenitors propose that muscle stem cells in a quiescent non-committed state exhibit enhanced engraftment potential (Hicks et al., 2018;Quarta et al., 2016;Montarras et al., 2005;Tierney et al., 2016).Our data demonstrate that upon activation and commitment dormant myogenic progenitors downregulate extracellular matrix proteins and upregulate expression of morphogens/receptors that make them susceptible to signals, like VEGFA, that communicates with vasculature during tissue reconstitution, or CD9, CD44, CD98 participating in cellular activation (Porpiglia et al., 2017;Figure 2-figure supplement 1B and D).Cell-cell communication analysis revealed that myogenic progenitors influence their own fate mainly with ECM-related signals .This finding is in line with the nature of in vivo fetal myogenic progenitors (Tierney et al., 2016) but also indicates that investigation at the myogenic progenitor level could provide new insights into ECM-related congenital muscular dystrophies.One example is the Ullrich congenital muscular dystrophy, where the ECM alters the muscle environment with progressive dystrophic changes, fibrosis, and evidence for increased apoptosis (Bönnemann, 2011).
CD82+ populations from our organoids engraft in the TA muscle of immunodeficient mice (Figure 4I  and J).It would be of interest for future studies to investigate whether increased engraftment can be achieved in 3D protocols (Faustino Martins et al., 2020;Shahriyari et al., 2022;ours) versus 2D patterned progenitor cells and to which degree this is attributed to high expression of extracellular matrix proteins.In particular, high Fibrillin1 expression on dormant non-dividing myogenic progenitors could potentially contribute to avoidance of fibrosis by myogenic progenitors through regulation of TGF-β signaling (Cohn et al., 2007).
Different phases of human fetal myogenesis have been modeled with 2D differentiation protocols (Shelton et al., 2014;Chal et al., 2015;Xi et al., 2017).Our 3D differentiation protocol does not go beyond these protocols when it comes to provide maturated physiologically responsive skeletal muscle cells, which we illustrate with the electrophysiological recording of organoid-derived cells of different origins (Figure 3K,L and M,.Structural distinctions like the posterior paraxial mesoderm on day 5 specified neural crest dermomytome on day 17, myogenic progenitor migration on day 23, and neural crest lineage arrest on day 35 (Figures 1C-H  and 2C) cannot be similarly observed in 2D protocols.In addition, our 3D organoid protocol provides myogenic progenitors in dormant and activated states for at least 14 wk in culture.We demonstrate that organoid culture sustains uncommitted MyoD1-negative, Pax7-positive myogenic progenitors as well as fibroadipogenic (PDGFRa+) progenitors, both resembling their fetal counterpart.This supply of muscle progenitors preserved in a quiescent state indicates translative potential of the approach.Future work will elucidate signaling transduction pathways during skeletal muscle organoid development to model and understand human myogenesis in more detail.LGMD2A-isogenic iPSC (Dorn et al., 2015;Mavrommatis et al., 2020;Panopoulos et al., 2017;Park et al., 2008) were cultured in TESR-E8 (StemCell Technologies) on Matrigel GFR (Corning)coated 6-well plates.Cells have been tested negative for mycoplasma.No experiments on human participants were performed.The use of reprogrammed human iPSC lines for research was performed after ethical approval from the ethics commission of the Ruhr-University Bochum, Medical Faculty (15-5401, August 2015).

Human skeletal muscle organoid differentiation protocol
Prior differentiation, undifferentiated human PSCs, 60-70% confluent, were enzymatically detached and dissociated into single cells using TrypLE Select (Thermo Fisher Scientific).Embryoid bodies formed via hanging drop approach with each droplet containing 3-4 × 10 3 human single PSCs in 20 μl were cultured hanging on TESR-E8 supplemented with polyvinyl alcohol (PVA) at 4 mg/ml (Sigma-Aldrich) and rock inhibitor (Y-27632) at 10 μM (StemCell Technologies) at the lid of Petri dishes.At the beginning of skeletal muscle organoid differentiation, embryoid bodies at the size of 250-300 μm were embedded into Matrigel and cultured in DMEM/F12 basal media (Thermo Fisher Scientific) supplemented with glutamine (Thermo Fisher Scientific), non-essential amino acids (Thermo Fisher Scientific), 100× ITS-G (Thermo Fisher Scientific), (Basal Media) 3 μM CHIR99021 (Sigma-Aldrich), and 0.5 μM LDN193189 (Sigma-Aldrich).On day 3, human recombinant basic fibroblast growth factor (bFGF) (PeproTech) at 10 ng/μl final concentration was added to the media.Subsequently, on day 5, the concentration of bFGF was reduced to 5 ng/μl and the media was further supplemented with 10 nM retinoic acid (Sigma-Aldrich).The differentiation media, on day 7, was supplemented only with human recombinant Sonic Hedgehog (hShh) (PeproTech) at 34 ng/μl, human recombinant WNT1A (PeproTech) at 20 ng/μl, and 0.5 μM LDN193189.On day 11, the cytokine composition of the media was changed to 10 ng/μl of bFGF and human recombinant hepatocyte growth factor (HGF) at 10 ng/ μl (PeproTech).From day 15 onward, the basal media was supplemented with ITS-X (Thermo Fisher Scientific) and human recombinant HGF at 10 ng/μl.The first three days of the differentiation the media were changed on a daily basis, from 3rd to 30th every second day, while from day 30 onward every 3rd day.The organoid approach was evaluated with six hiPSCs lines with independent genetic backgrounds, with more than five independent derivations per line, especially for the control line (CB CD34+) more than 20 derivation, obtaining always similar results.Per derivation and upon embryoid body Matrigel embedding, cultures exhibited high reproducibility.Upon migration and skeletal muscle formation, organoids are occupying whole Matrigel droplet and through generation of additional bulge reach sizes of 4-5 mM.Concomitantly, myogenic progenitors fall off the organoid and generate a sheath of myogenic organoids and muscle fibers at the surface of the culture plate.For all lines, functional myofibers and PAX7-positive myogenic populations could be evaluated.Myogenic populations from different lines exhibit high similarity (Pearson correlation, rho = 0.94-0,95, Figure 4-figure supplement 3G).

Cryosection immunochemistry
Organoids from different stages were fixed on 4% paraformaldehyde overnight at 4°C under shaking conditions, dehydrated (30% sucrose overnight incubation), and embedded in OCT freezing media.Cryosections were acquired on a Leica CM3050s Cryostat.For the immunostaining process, cryosections were rehydrated with PBS, followed by permeabilization once with 0.1% Tween-20 in PBS (rinsed 3× with PBS), and then with 0.1% Triton-X in PBS (rinsed 3× with PBS).Subsequently, the sections were blocked with 1% BSA/10% NGS in PBS for 1 hr at room temperature (RT).Primary antibody incubations were performed from 1 to 2 d at 4°C, with secondary antibody incubations for 2 hr at RT.

Oil O Red staining
For histological visualization of adipocytes within the organoids, Oil O Red Stain kit (Abcam, ab150678) was applied on frozen sections derived from PFA-fixated organoids following the manufacturer's recommended protocol.Organoid sections upon staining with Oil O Red were visualized with an Olympus BX61 upright microscope.

FACS EdU assay
At 15 wk post differentiation, organoids were incubated overnight with 5 μM EdU final concentration.The next day, organoids were dissociated into single cells by incubation at 37°C within papain solution for 1-2 hr, followed by incubation with TryplE Select for 10 min to ensure single-cell dissociation.Then, the dissociated cells were passed through a 70 μm cell culture strainer to remove any remaining aggregates.To detect EdU, the cells were processed with Click-iT EdU Alexa Fluor 488 Flow Cytometry Assay Kit (Invitrogen) according to the manufacturer's instructions and then analyzed using the 488 channel of a BD Biosciences FACSAria Fusion flow cytometer.
FACS isolation of ITGβ1 + /CXCR4 + myogenic cell population for RNA sequencing Organoids from Duchenne and Control iPSC lines and during 15th-16th week post differentiation were dissociated into single cells during incubation, with Papain solution dissociation upon gentle shaking (1-2 hr) was observed.To acquire singlets, the cells were filtered through 40 μm cell strainers and upon washing with 1% BSA solution prepared for surface antigen staining.For surface antigen staining, 20 min incubation with the Alexa Fluor 488 anti-human CD29 (BioLegend, clone TS2/16) and PE anti-human CD184[CXCR4] (BioLegend, clone 12G5) was applied together with the corresponding isotype controls: PE Mouse IgG2a, κ Isotype Ctrl antibody (BioLegend, clone MOPC-173) and 488 Mouse IgG1, κ Isotype Ctrl antibody (Invitrogen, clone P3.6.2.8.1).For removing residual antibodies, the cells were washed twice with 1% BSA staining solution and processed by BD Biosciences FACSAria Fusion flow cytometer.Briefly before FACS sorting to discriminate between dead and alive cells, DAPI was added to the samples and then DAPI -/CD29 + /CXCR4 + cell populations were collected into tubes containing RLT buffer supplemented with b-mercaptoethanol to avoid RNA degradation.The FACS gating strategy is further depicted in Figure 4-figure supplement 3A.

Bulk RNA sequencing RNA extraction
Total RNA was extracted from single organoids or cultured cells by using the RNAeasy Micro Kit (QIAGEN) according to the manufacturer's instructions.Subsequently, before library preparation, the RNA integrity was evaluated on an Agilent 2100 Bioanalyzer by using the RNA 6000 Pico kit (Agilent).cDNA library preparation: For 4 and 8 wk organoids, cDNA library was prepared by using the whole transcriptome Illumina TruSeq Stranded Total RNA Library Prep Kit Gold (Illumina), followed by evaluation on an Agilent 2100 Bioanalyzer by using the DNA 1000 kit.The resulting mRNA library was sequenced as 2 × 75 bp paired-end reads on a NextSeq 500 sequencer (Illumina).For 16 wk and ITGβ1 + /CXCR4 + sorted cells, cDNA library was prepared using the whole transcriptome Ovation Solo RNA seq Library Preparation Kit (TECAN, NuGEN), followed by evaluation on an Agilent 2100 Bioanalyzer by using the DNA 1000 Chip.The resulting mRNA library was sequenced as 1 × 150 bp single reads on a HiSeq 3000 sequencer (Illumina).

Bulk RNA-seq bioinformatic analysis
Sequenced reads were aligned to the human reference genome (hg38) with TopHat2 (version 2.1.1),and the aligned reads were used to quantify mRNA expression by using HTSeq-count (version 0.11.2).DESeq2 (Love et al., 2014) was used to identify differentially expressed genes (DEGs) across the samples.ITGβ1 + /CXCR4 + organoid-derived myogenic cell populations were compared to already available transcriptomic dataset of human fetal muscle progenitors (GSM2328841-2) (Hicks et al., 2018).

scRNA sequencing Sample and cDNA library preparation
Single cells were acquired upon incubation for 1 hr with solution containing papain and EDTA.Upon dissociation, the cell number and viability were estimated.Then, cells were resuspended in a solution containing 0.5% BSA in PBS to reach a concentration of 390 cells per μl.The cDNA library was prepared using the Chromium Single Cell 3′ Reagent Kits (v3): Single Cell 3′ Library & Gel Bead Kit v3 (PN-1000075), Single Cell B Chip Kit (PN-1000073), and i7 Multiplex Kit (PN-120262) (10x Genomics) according to the manufacturer's instructions.Then, the cDNA library was run on an Illumina HiSeq 3000 as 150 bp paired-end reads.

Single-cell RNA-seq bioinformatic analysis
Sequencing data were processed with UMI-tools (version 1.0.0),aligned to the human reference genome (hg38) with STAR (version 2.7.1a), and quantified with Subread featureCounts (version 1.6.4).Data normalization and further analysis were performed using Seurat (version 3.1.3, Stuart et al., 2019).For initial quality control of the extracted gene-cell matrices, cells were filtered with parameters low threshold = 500, high threshold = 6000 for number of genes per cell (nFeature_RNA), high threshold = 5 for percentage of mitochondrial genes ( percent.mito), and genes with parameter min.cells = 3. Filtered matrices were normalized by the LogNormalize method with scale factor = 10,000.Variable genes were found with parameters of selection.method = 'vst,' nfeatures = 2000, trimmed for the genes related to cell cycle (KEGG cell cycle, hsa04110) and then used for principal component analysis (PCA).Statistically significant principal components were determined by the JackStraw method and the first five principal components were used for nonlinear dimensional reduction (tSNE and UMAP) and clustering analysis with resolution = 0.2.Monocle3 (version 0.2.0,Cao et al., 2019) was used for pseudotime trajectory analysis.The data matrix of Seurat objects (assays[["RNA"]]@counts) was imported to the Monocle R package, then dimensionality reduction with the PCA method with parameters max_components = 2 was performed and then cluster_cells, learn_graph and order_cells functions were performed subsequently.Organoid-derived myogenic progenitors were compared to already available transcriptomic dataset of adult satellite cells (GSE130646) (Rubenstein et al., 2020).
For integrative analysis, we used the Seurat (version 4, Hao et al., 2021) package.For each dataset, cells with >6000 or <300 detected genes, as well as those with mitochondrial transcripts proportion >5-10%, were excluded.For finding anchors between control and DMD datasets using 5000 anchors and regressing out cell cycle genes, sequencing depth and stress-related genes were carried out before integration.Developmental score was calculated as described in Xi et al., 2020.Briefly, we used the ''AddMod-uleScore'' function to calculate the embryonic and adult score using a list of DEGs between adult and embryonic myogenic progenitor clusters.DEGs were selected from the Supplementary information in Xi et al., 2020 table mmc3.The developmental score was further calculated by subtracting embryonic from the adult score.Embryonic and fetal datasets were filtered using a low threshold = 500 and adult datasets were filtered using a low threshold = 250 for genes per cell.In addition, scaling the data (''S.Score,' 'G2M.Score,' 'Stress,' and 'total Count' to the vars.to.regress' argument) was used to regress out the effects of the cell cycle, dissociation-related stress, as well as cell size/sequencing depth to all datasets.Myogenic subpopulation of SMO and adult satellite cell datasets were selected by 'PAX7' expression.Embryonic and fetal (weeks 5-18) and adult satellite cell (years 7, 11, 34, and 42) scRNAseq data are from Xi et al., 2020 (GSE147457), adult satellite cell (year 25) scRNAseq data are from Rubenstein et al., 2020 (GSE130646).

Cell−cell communication analysis
To investigate cell−cell communications among activated, mitotic and dormant myogenic progenitors, fibroadipogenic and neural progenitors, and myofiber-related clusters from 12 wk organoids, the CellChat R package (Jin et al., 2021) was applied.

qPCR expression analysis
By pooling three organoids per sample, total RNA was extracted using the RNAeasy Mini Plus Kit (QIAGEN).For first-strand cDNA synthesis, the High Capacity RNA to cDNA Kit (Applied Biosystems) was applied using as template 2 μg of total RNA.For setting qPCR reactions, the GoTaq qPCR Master Mix (Promega) was used with input template 4 ng cDNA per reaction while the reaction was detected on a CFX 96 Real-Time PCR detection system (Bio-Rad).The relative quantification (ΔCT) method was applied for detecting changes in gene expression of pluripotent, neural tube, neural crest, posterior/anterior somitic, and dermomyotomal markers, between different time points along the differentiation.qPCR primers applied for each marker for evaluating organoid development are listed in Supplementary file 1a.

Diffusion map analysis
By pooling 3-5 organoids per sample, total RNA was extracted using the RNAeasy Mini Plus Kit (QIAGEN) and further processed as described under qPCR expression analysis.Normalized Ct values of the selected genes for each sample were used as input for generating eigenvector values using the destiny package (Angerer et al., 2016).Then, all samples were ordered by their Diffusion Component 1 at specific timepoints during early culture development (day 2-day 11).Ct values >30 were not considered for subsequent analysis.qPCR primers applied for each marker for diffusion map analysis are listed in Supplementary file 1b.

Second harmonic generation (SHG) imaging using multi-photon microscopy
A TriM Scope II multi-photon system from LaVision BioTec was used to visualize skeletal muscle fiber organization inside organoids and distinct sarcomeres.The microscope setup is a single-beam instrument with an upright Olympus BX51 WI microscope stand equipped with highly sensitive nondescanned detectors close to the objective lens.The TriM Scope II is fitted with a Coherent Scientific Chameleon Ultra II Ti:Sapphire laser (tuning range 680-1080 nm) and a Coherent Chameleon Compact OPO (automated wavelength extension from 1000 nm to 1600 nm).A 20× IR objective lens (Olympus XLUMPlanFl 20×/1.0W)with a working distance of 2.0 mm was used.Muscle fiber SHG signals were detected in forward direction using TiSa light at 850 nm, a 420/40 band pass filter, and a blue-sensitive photomultiplier (Hamamatsu H67080-01).3D images were acquired and processed with LaVision BioTec ImSpector Software.

Current measurement
Membrane currents were measured at ambient temperature (22-24°C) using standard whole-cell patch clamp software ISO2 (MFK, Niedernhausen, Germany).Cells were voltage-clamped at a holding potential of -90 mV, that is, negative to EnAChR, resulting in inward Na + currents.Every 10 s, voltage ramps (duration 500 ms) from -120 mV to +60 mV were applied to assess the stability of the recording conditions and generate I/V curves (membrane currents in response to depolarizing voltage ramps are shown as downward deflections).Signals were filtered (corner frequency, 1 kHz), digitally sampled at 1 kHz, and stored on a computer equipped with the hardware/software package ISO2 for voltage control, data acquisition, and data analysis.Rapid exposure to a solution containing acetylcholine was performed by means of a custom-made solenoid-operated flow system permitting a change of solution around an individual cell with a half-time of about 100 ms.For measurements, cells devoid of contact with neighboring cells were selected.Cells originated from organoids at week 8.

Fluorescence microscopy and imaging
To monitor changes in [Ca2 + ]i, skeletal muscle cells were transiently transfected with pcDNA3[Twitch-2B] (Addgene, 49531) (0.25 μg per 35 mm culture dish).Skeletal muscle cells were transfected using either poly-ethyleneimine (PEI) or Lipofectamine (Invitrogen) according to the manufacturer's instructions.Prior to experiments, cells were seeded on sterile, poly-L-lysine-coated glass coverslips and analyzed for 48 hr after transfections.All experiments were performed using single cells at ambient temperature.Fluorescence was recorded using an inverted microscope (Zeiss Axiovert 200, Carl Zeiss AG, Göttingen, Germany) equipped with a Zeiss oil immersion objective (100×/1.4), a Polychrome V illumination source, and a photodiode-based dual-emission photometry system suitable for CFP/ YFP-FRET (FEI Munich GmbH, Germany).For FRET measurements, single cells were excited at 435 nm wavelength with light pulses of variable duration (20-50 ms; frequency: 5 Hz) to minimize photobleaching.Corresponding emitted fluorescence from CFP (F480 or FCFP) or from YFP (F535 or FYFP) was acquired simultaneously, and FRET was defined as the ratio of FYFP/FCFP.Fluorescent signals were recorded and digitized using a commercial hardware/software package (EPC10 amplifier with an integrated D/A board and Patch-master software, HEKA, HEKA Elektronik, Germany).The individual FRET traces were normalized to the initial ratio value before agonist application (FRET/FRET0).

Transplantation experiments Flow cytometry
The organoid culture was dissociated using TrypLe (Thermo Fisher Scientific, 12563011) and filtered through 70-40 μm cell strainer.Single-cell suspension was stained with PE anti-human CD82 Antibody (BioLegend, 342103) in 2% BSA and 2 mM EDTA and sorted using FACS sorter (Beckman Coulter, MoFlo Astrios Zellsortierer).For gating, unstained single cells from the same organoid culture were used as a baseline control.Sorted cells were further processed for transplantation experiments.
The following dataset was generated: The following previously published datasets were used:

Figure supplement 1 .
Figure supplement 1. Lineage representation and organoid culture progression at early stages of differentiation protocol.

Figure supplement 2 .
Figure supplement 2. Pseudotime ordering of myogenic progenitor revealing distinct states and cell fate decisions.

Figure supplement 4 .
Figure supplement 4. Organoid-derived myogenic progenitors and correlation to adult human satellite cells.

Figure supplement 5 .
Figure supplement 5. Characterization of cell-cell communication network for all clusters at week 12 of human skeletal muscle organoid development.

Figure supplement 6 .
Figure supplement 6. Reproducibility of organoid culture at early and mature stages.

Figure supplement 7 .
Figure supplement 7. Comparison between 2D in vitro myogenic differentiation protocols and in vivo staging.