Self-organized fingering instabilities drive the emergence of tissue morphogenesis in digit organoids

Mechanochemical instabilities drive digit morphogenesis in organoids | bioRxiv /* */ /* */ <!-- <!-- /*! * yepnope1.5.4 * (c) WTFPL, GPLv2 */ (function(a,b,c){function d(a){return"[object Function]"==o.call(a)}function e(a){return"string"==typeof a}function f(){}function g(a){return!a||"loaded"==a||"complete"==a||"uninitialized"==a}function h(){var a=p.shift();q=1,a?a.t?m(function(){("c"==a.t?B.injectCss:B.injectJs)(a.s,0,a.a,a.x,a.e,1)},0):(a(),h()):q=0}function i(a,c,d,e,f,i,j){function k(b){if(!o&&g(l.readyState)&&(u.r=o=1,!q&&h(),l.onload=l.onreadystatechange=null,b)){"img"!=a&&m(function(){t.removeChild(l)},50);for(var d in y[c])y[c].hasOwnProperty(d)&&y[c][d].onload()}}var j=j||B.errorTimeout,l=b.createElement(a),o=0,r=0,u={t:d,s:c,e:f,a:i,x:j};1===y[c]&&(r=1,y[c]=[]),"object"==a?l.data=c:(l.src=c,l.type=a),l.width=l.height="0",l.onerror=l.onload=l.onreadystatechange=function(){k.call(this,r)},p.splice(e,0,u),"img"!=a&&(r||2===y[c]?(t.insertBefore(l,s?null:n),m(k,j)):y[c].push(l))}function j(a,b,c,d,f){return q=0,b=b||"j",e(a)?i("c"==b?v:u,a,b,this.i++,c,d,f):(p.splice(this.i++,0,a),1==p.length&&h()),this}function k(){var a=B;return a.loader={load:j,i:0},a}var l=b.documentElement,m=a.setTimeout,n=b.getElementsByTagName("script")[0],o={}.toString,p=[],q=0,r="MozAppearance"in l.style,s=r&&!!b.createRange().compareNode,t=s?l:n.parentNode,l=a.opera&&"[object Opera]"==o.call(a.opera),l=!!b.attachEvent&&!l,u=r?"object":l?"script":"img",v=l?"script":u,w=Array.isArray||function(a){return"[object Array]"==o.call(a)},x=[],y={},z={timeout:function(a,b){return b.length&&(a.timeout=b[0]),a}},A,B;B=function(a){function b(a){var a=a.split("!"),b=x.length,c=a.pop(),d=a.length,c={url:c,origUrl:c,prefixes:a},e,f,g;for(f=0;f<d;f++)g=a[f].split("="),(e=z[g.shift()])&&(c=e(c,g));for(f=0;f<b;f++)c=x[f](c);return c}function g(a,e,f,g,h){var i=b(a),j=i.autoCallback;i.url.split(".").pop().split("?").shift(),i.bypass||(e&&(e=d(e)?e:e[a]||e[g]||e[a.split("/").pop().split("?")[0]]),i.instead?i.instead(a,e,f,g,h):(y[i.url]?i.noexec=!0:y[i.url]=1,f.load(i.url,i.forceCSS||!i.forceJS&&"css"==i.url.split(".").pop().split("?").shift()?"c":c,i.noexec,i.attrs,i.timeout),(d(e)||d(j))&&f.load(function(){k(),e&&e(i.origUrl,h,g),j&&j(i.origUrl,h,g),y[i.url]=2})))}function h(a,b){function c(a,c){if(a){if(e(a))c||(j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}),g(a,j,b,0,h);else if(Object(a)===a)for(n in m=function(){var b=0,c;for(c in a)a.hasOwnProperty(c)&&b++;return b}(),a)a.hasOwnProperty(n)&&(!c&&!--m&&(d(j)?j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}:j[n]=function(a){return function(){var b=[].slice.call(arguments);a&&a.apply(this,b),l()}}(k[n])),g(a[n],j,b,n,h))}else!c&&l()}var h=!!a.test,i=a.load||a.both,j=a.callback||f,k=j,l=a.complete||f,m,n;c(h?a.yep:a.nope,!!i),i&&c(i)}var i,j,l=this.yepnope.loader;if(e(a))g(a,0,l,0);else if(w(a))for(i=0;i (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0];var j=d.createElement(s);var dl=l!='dataLayer'?'&l='+l:'';j.src='//www.googletagmanager.com/gtm.js?id='+i+dl;j.type='text/javascript';j.async=true;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-M677548'); Skip to main content Home About Submit ALERTS / RSS Search for this keyword Advanced Search New Results Mechanochemical instabilities drive digit morphogenesis in organoids View ORCID Profile Rio Tsutsumi , View ORCID Profile Antoine Diez , View ORCID Profile Steffen Plunder , Ryuichi Kimura , View ORCID Profile Oki Shinya , Kaori Takizawa , Rei Nakano , View ORCID Profile Haruhiko Akiyama , View ORCID Profile Ritsuko Takada , View ORCID Profile Shinji Takada , View ORCID Profile Marco Musy , View ORCID Profile James Sharpe , View ORCID Profile Mototsugu Eiraku doi: https://doi.org/10.1101/2025.08.31.673315 Rio Tsutsumi 1 Institute for Advanced Study of Human Biology (WPI-ASHBi), Kyoto University , Kyoto, 606-8303, Japan 2 Institute for Life and Medical Sciences (LiMe), Kyoto University , Kyoto, 606-8507, Japan Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Rio Tsutsumi For correspondence: tsutsumi.rio.8x{at}kyoto-u.ac.jp eiraku.mototsugu.2u{at}kyoto-u.ac.jp Antoine Diez 1 Institute for Advanced Study of Human Biology (WPI-ASHBi), Kyoto University , Kyoto, 606-8303, Japan 3 Mathematical Application Research Team, Division of Applied Mathematical Science, RIKEN Center for Interdisciplinary Theoretical and Mathematical Sciences (iTHEMS) , RIKEN iTHEMS Wako, Saitama 351-0198, Japan Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Antoine Diez Steffen Plunder 1 Institute for Advanced Study of Human Biology (WPI-ASHBi), Kyoto University , Kyoto, 606-8303, Japan Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Steffen Plunder Ryuichi Kimura 4 Department of Drug Discovery Medicine, Graduate School of Medicine, Kyoto University , Kyoto, 606-8507, Japan 5 Institute of Resource Development and Analysis, Kumamoto University , Kumamoto, 860-0811, Japan Find this author on Google Scholar Find this author on PubMed Search for this author on this site Oki Shinya 4 Department of Drug Discovery Medicine, Graduate School of Medicine, Kyoto University , Kyoto, 606-8507, Japan 5 Institute of Resource Development and Analysis, Kumamoto University , Kumamoto, 860-0811, Japan Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Oki Shinya Kaori Takizawa 2 Institute for Life and Medical Sciences (LiMe), Kyoto University , Kyoto, 606-8507, Japan Find this author on Google Scholar Find this author on PubMed Search for this author on this site Rei Nakano 2 Institute for Life and Medical Sciences (LiMe), Kyoto University , Kyoto, 606-8507, Japan 6 Department of Polymer Chemistry, Graduate School of Engineering, Kyoto University , Kyoto, 615-8510, Japan Find this author on Google Scholar Find this author on PubMed Search for this author on this site Haruhiko Akiyama 7 Department of Orthopaedic Surgery, Graduate School of Medicine, Gifu University , Gifu, 501-1194, Japan 8 Center for One Medicine Innovative Translational Research (COMIT), Institute for Advanced Study, Gifu University , Gifu, 501-1193, Japan Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Haruhiko Akiyama Ritsuko Takada 9 National Institute for Basic Biology, National Institutes of Natural Sciences , Aichi, 444-8787, Japan 10 Exploratory Research Center on Life and Living Systems (ExCELLS), National Institutes of Natural Sciences , Aichi, 444-8787, Japan Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Ritsuko Takada Shinji Takada 9 National Institute for Basic Biology, National Institutes of Natural Sciences , Aichi, 444-8787, Japan 10 Exploratory Research Center on Life and Living Systems (ExCELLS), National Institutes of Natural Sciences , Aichi, 444-8787, Japan Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Shinji Takada Marco Musy 11 European Molecular Biology Laboratory (EMBL-Barcelona) , 08003, Barcelona, Spain 12 Elisava, Barcelona School of Design and Engineering, University of Vic-Central University of Catalonia (UVIC-UCC) , Barcelona, Spain Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Marco Musy James Sharpe 11 European Molecular Biology Laboratory (EMBL-Barcelona) , 08003, Barcelona, Spain 13 Barcelona Collaboratorium for Modelling and Predictive Biology , 08005 Barcelona, Spain 14 Institució Catalana de Recerca i Estudis Avançats (ICREA) , 08010 Barcelona, Spain Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for James Sharpe Mototsugu Eiraku 2 Institute for Life and Medical Sciences (LiMe), Kyoto University , Kyoto, 606-8507, Japan 6 Department of Polymer Chemistry, Graduate School of Engineering, Kyoto University , Kyoto, 615-8510, Japan Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Mototsugu Eiraku For correspondence: tsutsumi.rio.8x{at}kyoto-u.ac.jp eiraku.mototsugu.2u{at}kyoto-u.ac.jp Abstract Full Text Info/History Metrics Supplementary material Preview PDF Abstract The emergence of complex anatomical structures —such as the hands— from unstructured tissues remains a fundamental question in developmental biology. Turing-type reaction-diffusion models have provided a molecular explanation for the periodic pre-patterning of digits; however, the physical principles driving 3D morphogenesis remain incompletely understood. To identify the biophysical design principles leading to digit formation, we develop a limb-mesenchymal organoid system that spontaneously forms elongated, digit-like protrusions. Iterations between experiments and agent-based models at the cellular level identify sufficient microscopic mechanisms leading to morphogenesis of digit-like structures: symmetry-breaking and the elongation of digits result from a combination of differential cell adhesion and morphogen-induced chemotaxis and convergent-extension. Lastly, to describe tissue-scale deformations, we perform a coarse-graining analysis of the agent-based model and derive a continuum model that reveals a structural analogy to Cahn-Hilliard-type equations. These equations are typically used to describe fluid phase separation and so-called « fingering instabilities » in fluid physics. Here, we show that they also accurately describe organoid morphogenesis. These findings suggest that “finger” formation is driven by a mechanical fingering instability acting in concert with chemical patterning, shedding a new light on vertebrate limb morphogenesis. Main The morphology of anatomical structures has been evolved to execute specialized functions with remarkable efficiency. A prime example is the distal end of the vertebrate limb, which consists of a modular design of repetitive elements that form the digits, enabling complex tasks such as goal-directed manipulation. Although these structures vary in shape, size, and number across species, their core design principles —governed by molecular and cellular behaviors during development— are fundamentally conserved 1 . Despite our detailed knowledge of the molecular and cellular behaviors involved, a fundamental question remains: how do microscopic interactions translate into the macroscopic shape of tissues? 2 This underscores the emergent nature of morphogenesis, where tissue-level structures cannot be predicted solely from a list of individual components 3 . In physics, similar problems have been addressed through statistical physics frameworks that bridge microscopic interaction rules and macroscopic laws. Whether tissue morphogenesis can be understood in an analogous way is an open and challenging question that requires resolving the complexity of experimental systems into analytically tractable multiscale mathematical models. In this context, a key concept is the idea of instability, which postulates that complex feedback loops can amplify small microscopic perturbations to produce macroscopic structures. One of the most influential instabilities linking molecular behaviors and tissue-scale pattern formation in biology is the Turing instability, driven by biochemical reaction–diffusion systems 4 . It effectively explains early patterning of digit and interdigit primordial tissues (the latter subsequently undergoes apoptosis) in the mouse limb bud (E10.5), where patterns emerge prior to significant cell movement 5 , 6 . However, while sufficient for initial molecular symmetry-breaking, this framework fails to account for the tissue-scale morphogenetic movements that physically shape digits along the distal–proximal axis by E12.5. The morphogenesis at later stages should depend on mechanochemical cell–cell interactions 7 , yet whether macroscopic morphogenesis can emerge from microscopic cellular behaviors through a mechanochemical instability remains an open question. To address this, we seek to identify the minimal biomechanical and biophysical "design principles" essential for digit formation. Multicellular in vitro systems, particularly organoids, provide a powerful platform to study such questions, as they can self-organize and give rise to multicellular patterns and structures under defined controllable conditions, thereby bypassing much of the complexity of in vivo development 8 – 11 . Here, we establish an in vitro organoid model of digit morphogenesis and combine it with agent-based and continuum mathematical models to identify the minimal design principles sufficient for digit formation. We propose the existence of a mechanochemical instability mechanism capable of driving tissue-level morphogenesis as an inevitable macroscopic consequence of microscopic cellular interaction rules, thus providing a framework for how molecular and biophysical components collectively give rise to the emergence of anatomical shapes. Fgf8b and Wnt3a signalling are sufficient to induce digit-like structures in limb mesenchymal organoids To identify the fundamental design principles required for digit morphogenesis of the hand, we first dissociated and cultured 9,000 limb mesenchyme cells derived from Sox9-GFP mouse embryos (to visualize cartilage differentiation), which formed spherical cellular aggregates in a three-dimensional (3D) culture conditions ( Fig. 1A ). Aggregates cultured under serum-only conditions (negative control: NC) largely maintained their spherical shape, but displayed non-uniform Sox9 expression ( Fig. 1B , left panel, Supplementary movie 1). As Fgf8b and Wnt3a are known to be secreted from the apical ectodermal ridge (AER) 12 and the overlying ectoderm, respectively, and act synergistically to induce distal identity and growth of limb mesenchyme morphogenesis 13 – 15 , we next investigated the effect of Fgf8/Wnt3 (FW) on cultured mesenchymal cellular aggregates. We observed that FW cultured aggregates not only went beyond non-uniform Sox9 distribution, but also broke their shape-symmetry growing into stick-shaped cartilage elongations, despite the fact that we supplemented FW uniformly in the medium ( Fig. 1B middle panel). Interestingly, when cultured with a larger cellular aggregate (27,000 cells), we observed the formation of multiple bifurcated cartilage rather than an isotropic scale up ( Fig. 1B right panels, 1C). This indicates that FW is sufficient to initiate not only shape-symmetry breaking, but also the physical elongation of cartilage without any directional cues in the medium ( Fig. 1D ). Download figure Open in new tab Figure 1. Formation of digit-like structures in organoid cultures. (A) Schematic: limb mesenchyme from Sox9–KI–EGFP E11.5 mouse embryos is aggregated and cultured. (B) Endpoint morphology at 193 h in NC or FW (Fgf8+Wnt3a) conditions, initiated with 9,000 or 27,000 cells; Sox9-GFP marks cartilage (green) on bright-field (scale bars, 1 mm). (C) The number of protrusions depending on the starting number of cells. Stacked bar charts show the percentage of colonies with a given number of protrusions across starting-cell groups; segment labels are raw counts and x-axis ticks show total n per group. Association between starting cell number and protrusion-category distribution was tested using a batch-stratified Cochran–Mantel–Haenszel (CMH) test, restricted to 2 batches containing ≥2 starting-cell groups: CMH M 2 = 28.07, df = 10, p = 0.00176 (two-sided, α = 0.05). Percentages in bars are descriptive. (D) Time-lapse in FW from 9,000 cells for 193 h, showing protrusion emergence and elongation (scale bar, 1 mm). (E) Distal–proximal lineage labeling. The cells from the distal half of autopod are labeled H2B-Venus (green) and proximal cells H2B-tdTomato (magenta) before aggregation, and cultured in an organoid. (F) Labeled organoids at 76 h in NC and FW, started from 9,000 or 27,000 cells, overlaid with bright-field (scale bar, 500 µm). Tips were enriched in Distal cells. To determine the mechanistic basis of symmetry-breaking in our aggregates, we explored a hypothesis that it is based on the cell sorting between two cell populations —the cells originated from distal (digit) and proximal (palm) parts of the autopod. We dissected and labelled distal and proximal cells with H2BVenus and H2BtdTomato, respectively and cultured the aggregates of the two cells under NC and FW conditions ( Fig. 1E ). Although distal and proximal cells formed an intermixed cellular aggregate under NC condition, they underwent sorting in FW conditions within 24 hrs, with the distal cell cluster subsequently forming an elongating tip region ( Fig. 1F , left and middle panels, Supplementary movie 2). This suggests that shape asymmetry is dependent on cell sorting . Interestingly, when more cells were used as starting material (27,000 cells), distal cells tended to cluster at multiple sites, indicating the bifurcation of cartilage ng ( Fig. 1F , right panel). Next, we performed a high-depth transcriptome analysis for a locally defined area via photo-isolation chemistry (PIC) 16 , to identify the genes highly upregulated at the tip of the organoids. Combined with our scRNA-seq analysis the results showed not only that the tips of organoids corresponded well with the genes expressed in the distal digits, including Hoxd13, Msx1, Lpatternihx2, and Sox11 17 – 20 ( Fig. S1A ), but also suggested that FW maintains distal and proximal identities, while NC drives both distal and proximal cells towards proximal identity and into an “inactive” state, aligning with the “Progress Zone model” of limb development 21 (Supplementary_scRNA-seq, Fig. S1 and S2 ). Therefore, based on our cell lineage and expression profile data, this suggests that FW is sufficient to trigger the emergence of a digit-like morphological structure in our 3D limb-mesenchymal organoid system, thereby mimicking faithfully the digit formation in vivo . Download figure Open in new tab Figure S1. Distal–Proximal gene programs in digit organoids (FW). (A) Volcano plot comparing tip (Distal) vs bottom (Proximal) PIC profiles; x = log₂FC, y = −log₁₀(p). Genes called Tip (red) if log₂FC > 0.6 & p < 0.05, Bottom (blue) if log₂FC < −0.6 & p < 0.05; others black. Vertical lines ±0.6 (∼1.5-fold); horizontal line p = 0.05; significant genes labeled. (B) UMAP of scRNA-seq: Day 1 (red), Day 2 (blue). (C) Marker table for surface ectoderm, AER, limb mesenchyme: % cells with ≥1 UMI; summary row shows % cells with ≥2 markers/tissue. (D) DotPlot of these markers across clusters (Day1/Day2 Distal/Proximal; muscles): dot size = %≥1 UMI; color = scaled mean. (E–F) UMAP FeaturePlots of Distal_Score (average of top 20 tip-enriched genes from A; magenta→green) and Proximal_Score (top 20 bottom-enriched; green→magenta). (G) UMAP annotated by Day and Distal/Proximal using B, E–F. (H) RNA-velocity streamlines (scVelo). (I) Slingshot principal curves, rooted as in H, with endpoints Day2_FW_Distal_2 (lineage 1: Distal) and Day2_FW_Proximal (lineage 2: Proximal). (J) For selected genes (Wnt5a, Fgf10, Inhba, Bmp2/4, Hoxd13): top, UMAP colored by log(expression+1) with Slingshot curves; bottom, log(expression+1) vs Slingshot pseudotime per lineage (Distal: purple, Proximal: yellow) with tradeSeq NB-GAM smoothers. Download figure Open in new tab Figure S2. FW condition maintains Distal identity and an “active” state. (A) UMAP of digit organoids (FW/NC; Day 1/Day 2) together with E11.5 Distal/Proximal cells (Day 0), annotated by DNA–lipid barcoding and Fig. S2G . (B–C) FeaturePlots of per-cell module scores on the UMAP. The FW signature is the average expression of the top 50 genes up-regulated in Day1_FW vs Day2_NC (Wilcoxon FindMarkers, padj < 0.05); the NC signature is the corresponding NC-up set. Scores are shown per cell with truncated scales (10th–90th percentiles): (B) FW signature (green→magenta), (C) NC signature (magenta→green). (D) Violin/box plots of the FW signature across binned clusters (Day1-/Day2-FW-Distal, Day1-/Day2-FW-Proximal, Day1-/Day2-NC); white dots mark medians. Pairwise two-sided Wilcoxon rank-sum tests; Benjamini–Hochberg–adjusted p-values are shown above brackets. x-axis labels include per-bin sample sizes (n); values are module-score units. (E) Schematic: FW sustains Distal identity and an “active” program. (F–G) Metacell scatter 64 : x = FW module score, y = NC module score. (F) Colored by metacell class (P1 naïve, P3 autopodial, P2 proximal progenitors; CT0–CT2 connective tissue/tendon; C1–C6 chondrocyte states; Attachment). (G) Same scatter colored by developmental time (E10.5–E14.5; majority vote per metacell). To analyze the contribution of cell sorting to the patterning of Sox9 expression in quantitative manners, we also cultured H2BtdTomato-labelled distal and H2BtagBFP-labelled proximal cells, from the Sox9-GFP mouse embryos on two-dimensional (2D) micromass culture ( Fig. S3A ). Similar to our findings in limb-mesenchymal organoids, we observed that FW conditions promoted cell sorting between the distal and proximal cells in comparison to NC conditions, ( Fig. S3B , Supplementary movie 3). Interestingly, when spots of Sox9 expression emerged, they typically appeared to be associated with either sorted clusters of distal cells or proximal cells. This is consistent with limb development in vivo , as the molecular pre-patterning of Sox9 corresponds not only to distal digital elements, but also the more proximal metacarpal and carpal skeletal elements ( Fig. S3B ). Strikingly, although Sox9 spots could form in distal clusters or proximal clusters, they appear less likely to emerge in regions with strong intermixing of proximal and distal (or on the border between the cluster), indicating that the reaction-diffusion mechanism of molecular patterning operates strongest between cells of the same proximo-distal identity. Download figure Open in new tab Figure S3. Variations of digit organoid cultures. (A) Schematics: Distal–proximal lineage labeling. The cells from the distal and proximal half of autopod of Sox9–KI–EGFP (white) E11.5 mouse embryos are labeled before aggregation, and cultured as micromass culture. (B) Micromass culture of limb mesenchyme in NC and FW at 72 h; Sox9-GFP in white; Distal and Proximal cells labeled with H2BtdTomato (green) and H2BtagBFP (magenta). (C) Dot plots of Sox9-GFP intensity over time in Distal (left) and Proximal (right) cells; points colored by the number of Distal neighbors among the 10 nearest. (D) Violin/box plots at 72 h comparing neighbor groups (≥9 vs ≤1 Distal neighbors) within each cell type using two-tailed unpaired Student’s t-tests; significance shown as asterisks ( p <0.05, * p <0.01, ** p <0.001). Percentages per group are shown below; sample sizes ( n ) and means are annotated. (E) Average Nearest Neighbor Score (NNS, k=10) over time in NC and FW, defined as the proportion of same-type neighbors normalized by cell-type frequency; lines show means, ribbons show bootstrapped 95% CIs. To quantify the level of cell sorting, we computed the nearest neighbor score (NNS) by counting the number of the same type of cells (distal or proximal) and normalizing it such that the score would be one under the null hypothesis (random distribution), while the score would be higher when cell sorting occurs (See methods) ( Fig. S3C ). This confirmed that the proportion of distal cells with at least nine distal neighbours (out of 10 neighbouring cells) increased over time, and those cells exhibited higher Sox9-GFP intensities. This effect was weaker for proximal cells, but nevertheless consistent: over time the more clustered proximal cells increase their expression of Sox9 ( Fig. S3D, E ). These results suggest that under FW conditions, the cell sorting between distal and proximal cells has a strong influence on the spatial patterning of cartilage ( Fig. S3D,E , S4A ). In contrast, when only distal cells were cultured in NC condition, labyrinth-like patterns were created, and the intervals between these nodules were further extended under FW conditions, which has been attributed to Fgf wavelength modification ( Fig. S4B, C ) 5 , 6 , 22 . Download figure Open in new tab Figure S4. 2D micromass culture in various cell types and culture conditions. (A) Distal + Proximal cells in FW, (B) Distal cells in FW, (C) Distal cells in NC at 72 h. DAPI (blue), Distal (green), Proximal (magenta), Sox9–GFP (white). Scale bar, 500 µm. Cell adhesion and chemotaxis drives symmetry-breaking but not digit-like structure elongation Considering that dynamic cell movements are required for symmetry-breaking and elongation observed in our 3D limb-mesenchymal organoids, we next built an Agent (cell)-Based Model (ABM). Cell sorting is typically achieved through differential cell adhesion. Indeed, limb-mesenchymal cells at different positions and ages have different cell adhesiveness 23 . For example, EphA7 is known to regulate the cell sorting between distal and proximal autopod mesenchyme at E11.5 mouse embryos 24 . Cell adhesiveness can be parametrised as the contact angle between cell doublets, whereby cells with stronger adhesion have larger angles and vice versa 25 ( Fig. 2A , left schematic). Thus, we cultured aggregates of distal and proximal cells separately under NC and FW conditions, dissociated on them on day 1 and day 2 to make doublets, and measured the angles between cells ( Fig. S5A ). We found that cell adhesions between distal-distal cell interactions were stronger than distal-proximal and proximal-proximal on day 1, and that the distal-proximal adhesion became equally strong as distal-distal on day 2 when cultured in FW, while no differences in cell adhesions were observed when cultured in NC ( Fig. 2A , Fig. S5B ). Download figure Open in new tab Figure S5. Details on the mechanisms of the morphological symmetry breaking. (A) Doublet assay workflow. Distal and Proximal cell aggregates were cultured in FW, dissociated on Day 1 or Day 2, and membrane-labeled with MemGlow dyes (488 or 560). Doublets were then assembled in predefined pairings: Distal_488–Distal_560 (DD), Distal_560–Proximal_488 (DP), and Proximal_488–Proximal_560 (PP). (B) Beeswarm+box plots of Angle by doublet class (DD, DP, PP) for Day1/Day2 NC. Dots = doublets; shape = batch. Stats: one-way ANOVA per experiment; if significant (α=0.05), Tukey HSD shown. (C) Agent-based model with variable parameters of differential adhesion (DD > DP > PP) and a traction bias toward Fgf8b (outer→inner). Axes: α, adhesion contrast; η, traction bias. Dotted box: adhesion-only. Solid box: adhesion + traction with an inset showing outline of a protrusion (dotted). (D) Live imaging of assembly assay at 96h: a proximal signaling aggregate (Wnt5a-P2A-tagBFP or Inhba-P2A-tagBFP) assembled with a Distal+Proximal aggregate. (E) Signal distribution at 96 h. Intensities for tagBFP, H2B-Venus (Distal), and H2B-tdTomato (Proximal) were integrated along a line from the tagBFP region; intensity-weighted medians defined peaks. Distances from the tagBFP peak to Distal vs Proximal peaks were compared by paired one-sided Wilcoxon (hypothesis: Distal closer). Plots show mean normalized profiles ±95% CI with median markers. Download figure Open in new tab Figure 2. Differential adhesion plus Fgf8b-directed traction bias explains symmetry breaking. (A) Beeswarm+box plots of Angle by doublet class (DD, DP, PP) for Day1/Day2 FW. Dots = doublets; shape = batch. Stats: one-way ANOVA per experiment; if significant (α=0.05), Tukey HSD shown. (B) Agent-based model implementing differential adhesion only (DD>DP>PP) (C) Day-1 FW immunostaining: dpERK and 7Tcf-EGFP (red→blue). (D) Rule map: angle to the Fgf8b gradient (x-axis) vs traction strength (y-axis) . DD (green, strong) and PP (magenta, weak) are constant; DP (yellow) is angle-dependent, biasing toward Fgf8b gradient. Dotted box: DP vectors (blue arrows) drive Distal motion (black arrowhead); Solid box: predicted Distal motion toward Fgf8b. (E) Agent-based model implementing differential adhesion + traction bias toward the Fgf8b gradient. Right: Fgf8b field. (F) Assembly assay: a proximal signaling aggregate (tagBFP or Fgf8b-P2A-tagBFP) assembled with a Distal+Proximal aggregate. (G) Live imaging at 96 h. (H) Signal distribution at 96 h. Intensities for tagBFP, H2B-Venus (distal), and H2B-tdTomato (proximal) were integrated along a line from the tagBFP region; intensity-weighted medians defined peaks. Distances from the tagBFP peak to Distal vs Proximal peaks were compared by paired one-sided Wilcoxon (hypothesis: Distal closer). Plots show mean normalized profiles ±95% CI with median markers. Based on these results, we constructed a minimal representative ABM that takes into account the key biophysical mechanisms leading to the emergence of digits from a mixed aggregate of distal and proximal cells. By incorporating differential cell adhesion, whereby distal-distal adhesion is stronger than the distal-proximal or proximal-proximal, we found that distal cells only cluster at the centre of the organoid ( Fig. 2B , S5C, Supplementary movie 4). This indicates mechanistically that differential cell adhesion alone is insufficient for distal cells to cluster on the surface of our organoids and that an additional mechanism is also likely to be involved. As exposure to FW in the culture medium may lead to an inherent bias in signalling cues received between the surface and centre of the organoid, we investigated the activation of the downstream signalling pathways in response to FW within the organoids. We performed immunostaining against phosphorylated-ERK and electroporated 7Tcf-GFP, to visualize the activation of Fgf signalling in response to Fgf8b and canonical Wnt signalling in response to Wnt3a, respectively. Indeed, in both cases, we observed that cells at the surface of the organoid had stronger signal activation than at the centre ( Fig. 2C ). Given that a gradient of FW is inherently established within the organoid towards its surface, we next tested whether chemotaxis of distal cells towards the FW signals may contribute to symmetry-breaking and digit elongation. Here, we modelled cell motion solely with reciprocal pulling/pushing mechanisms between pairs of cells that preserves the symmetry of Newton’s third law 26 . For directed cell migration towards the gradient in the ABM, in addition to the stronger adhesion strength of distal cells, we also implemented that distal cells have a directional traction bias. This was defined such that when distal cells are adjacent to proximal cells, they pull the proximal cells and propel themselves strongly towards the chemotactic gradient (Fgf8b and/or Wnt3a, in this case) ( Fig. 2D ). By using such parameters, the ABM was able to recapitulate the morphological symmetry-breaking in limb-mesenchymal organoids, but not elongation beyond spherical clusters of distal cells located at the surface ( Fig. 2E , dotted line, S5C, Supplementary movie 5). To further investigate the chemotactic behaviour of distal cells towards FW, we generated limb-mesenchymal organoids that included small aggregates of proximal cells which overexpressed either Fgf8b-T2A-tagBFP or tagBFP, as a control ( Fig. 2F ). As expected, after culturing under NC conditions, distal cells moved to the close proximity of the Fgf8b-T2A-tagBFP overexpressing cells, but not to tagBFP overexpressing cells ( Fig. 2G, H , Supplementary movie 6). Importantly, chemotactic activity was not observed when known morphogens expressed in limb mesenchyme, such as Wnt5a or Inhba (encoding activin), were overexpressed in these small aggregates ( Fig. S5D, E ). These results showed that distal cells showed chemotaxis towards at least Fgf8b gradient. Wnt5a regulates sustained digit elongation via convergent extension Although our ABM successfully reconstructed both the shape symmetry-breaking events, it was unable to recapitulate the digit elongation process beyond spherical distal clusters. Since both Fgf8 and Wnt3 have been shown to promote cell proliferation in a context-dependent manner 13 , we investigated the role of FW in inducing cell proliferation in our organoid system and whether this contributes to symmetry-breaking and digit elongation. We generated limb-mesenchymal organoids and performed EdU-pulse labelling to determine how NC and FW conditions affect distal and proximal cell proliferation, after two hours of chasing on Day 1 and 2. Although FW conditions promoted a higher proliferation rate in comparison to NC conditions, irrespective of the stage, the ratio of distal to proximal cells remained relatively unchanged ( Fig. 3A ). Furthermore, when cell proliferation was inhibited using Trichostatin A (TSA), tissue elongation was still observed, albeit as smaller aggregates ( Fig. S6A ). These results indicate that an alternative cellular mechanism is required for the sustained elongation of the distal cell clusters. Download figure Open in new tab Figure S6. Tissue elongation is proliferation-independent, and Wnt5a is enriched in Distal cells. (A) Inhibitor treatment for cell proliferation. Experimental schedule (left): organoids were cultured in FW (Fgf8b/Wnt3a) from Day0–Day2; TSA (HDAC inhibitor) was added Day1–Day2; BrdU was pulsed Day1-Day2. Bright-field images at Day1 and Day2 are shown for 0 nM and 20 nM TSA, with BrdU immunostaining (magenta) on Day2 (right). scale bars, 300 µm. (B) Relative Wnt5a expression (qPCR) in 2D micromass, Day 1. Box–dot plots display Day1_FW_Distal, Day1_FW_Proximal, Day1_NC_Distal, and Day1_NC_Proximal. Each dot is one biological replicate; the y-axis shows expression normalized to Gapdh. Statistics: one-way ANOVA across the four groups with Tukey’s HSD for pairwise comparisons. Planned Distal vs Proximal contrasts within FW and within NC were additionally tested by two-sided t-tests; the unadjusted p values are printed on the plot. Download figure Open in new tab Figure 3. Wnt5a expression in Distal cells and convergent extension. (A) EdU incorporation (% EdU⁺) in Distal (green) vs Proximal (magenta) cells on Day 1/Day 2 under NC and FW. Bars show means; dots are individual samples (shape = Batch). Two-sided unpaired t-tests (Distal vs Proximal); unadjusted p values shown. (B) Day-2 FW organoid immunostaining of Wnt5a (green) with DAPI (blue). (C–E) Convergent-extension (CE) analysis. (D) Workflow: tracked cell positions are smoothed into a velocity field and integrated to estimate tissue deformation; from this, the primary axis of extension and its strength are extracted and visualized as deformation of a reference spheres into ellipses (See Supplementary document). (E) Representative CE analysis from experimental data. (F) CE vectors projected along the tissue-elongation axis are plotted versus positions on the tissue-elongation axis (red) and overlaid with Distal (green) and Proximal (magenta) cell-density profiles, highlighting regions where CE occurs. (F) Overexpression screen in Proximal cells. (Left) Schematics and example of the organoids. (Right) Box–jitter plots show organoid aspect ratio by gene (x); each dot is one aggregate (shape = Batch). Control = tagBFP-only; constructs express Bmp2, Bmp4, Fgf10, Hoxd13, Inhba, Wnt5a. Statistics: planned two-sided t-tests vs Control; unadjusted p values are printed above brackets. Tissue elongation can occur without exerting force on an external substrate via push-pull mechanisms while maintaining symmetrical cell-cell interactions, by undergoing the process of convergent extension. This can be achieved through the establishment of a morphogen gradient that allows cells to pull each other orthogonally along the gradient or axis of elongation 26 . Previous studies have shown that Wnt5a triggers planar cell polarity, which is known to be necessary for convergent extension in the digit 27 – 29 . Indeed, we observed that distal cells expressed higher levels of Wnt5a in 2D micromass culture and at the distal tips in 3D limb mesenchymal organoids ( Fig. S6B , Fig. 3B ). This suggests that the accumulation of morphogens such as Wnt5a in distal clusters could potentially trigger convergent extension in our organoid system. To determine whether convergent extension contributes to digit elongation in our limb-mesenchymal organoids, we performed two-photon live-imaging microscope 11 to create four-dimensional (4D) images and tracked the cell trajectories using TrackMate 30 and ELEPHANT 31 . Additionally, we developed a novel imaging data analysis pipeline to detect and quantify tissue-level deformations. First, we interpolated the velocity vector fields obtained from the tracks using a Gaussian kernel smoothing method. From the smoothed vector fields, we computed a deformation map of the organoid based on a notion originating from continuum mechanics, called Cauchy-Green strain tensor. This tensor measures how a reference sphere in the initial frame of the experiment would be deformed by the vector field. This deformation is represented by ellipses, whereby the longer axes indicate the directions of tissue elongation ( Fig. 3C-E ). We then plotted the primary axes of deformation of these ellipses and projected them onto the elongation axis of the digit organoid, which essentially measures the extent of convergent extension (see Supplementary_CauchyGreenTensorEstimation). Our analysis showed indeed that convergent extension is maximal where the Wnt5a gradient is the steepest, at the distal-proximal interface ( Fig. 3E ). As gradient formation is required for convergence extension, to perturb this distal-proximal gradient, we electroporated proximal cells with the distal morphogen gene construct CAG-tagBFP-2A-Wnt5a as well as other potential distal morphogen genes, including Bmp2, Bmp4, Fgf10, and Inhba (encoding Activin). Additionally, since Hoxd13 is the upstream of the distal identities, we used it as a positive control 32 , 33 . We found that organoids expressing Hoxd13, Inhba, and Wnt5a showed a significant reduction in the aspect ratios (the longest axis vs orthogonal axis; Fig. 3F , Supplementary movie 7). In vivo , the distal-proximal gradient of Wnt5a controls convergent-extension and promotes the formation of the expression domain of activin 29 , 34 . Together, these findings indicate that a distal-proximal gradient (eg. Wnt5a) is necessary to establish convergent extension for sustained digit elongation. The minimal sufficient mechanisms for digit morphogenesis in limb-mesenchymal organoids Based on our data, we added another traction bias in our ABM defined such that two adjacent distal cells pull each other more strongly in the direction orthogonal to the distal-proximal elongation gradient (Wnt5a in this case). In combination with a standard (symmetrical) repulsion force between pairs of cells, this mechanism then allows for a tissue-level convergent-extension motion ( Fig. 4A ). The ABM incorporating differential adhesions, distal expression of Wnt5a, chemotaxis, and convergent-extension successfully recapitulated the symmetry-breaking and elongation, showing that these mechanisms are sufficient to explain digit morphogenesis in limb-mesenchymal organoids ( Fig. 4B , see dotted line, Supplementary movie 8, Fig. S7A ). Download figure Open in new tab Figure S7. Tissue elongation is dependent on traction bias orthogonal to Wnt5a gradient. (A) ABM combining differential adhesion (DD>DP>PP), DP bias toward Fgf8b, Distal Wnt5a production, and DD bias orthogonal to Wnt5a. Parameters: traction bias ( η) vs adhesion differences ( α ) solid box marks protrusion shape (dotted line in inset). (B) CE index in an ABM with differential adhesion and Fgf8b-directed traction bias. The solid red curve shows the CE index projected onto the tissue-elongation (AP) axis; the gray dashed curve shows the orthogonal component. Shaded green and magenta bands are kernel-density estimates of Distal and Proximal cells, respectively (four time points). The AP-aligned CE peak is modest and spatially coincident with Distal cluster, unlike the larger peaks at the interface between Distal and Proximal clusters in Figs. 3F and 4D where Wnt5a-induced/orthogonal biases are included. Download figure Open in new tab Figure 4. Sufficient conditions to explain tissue deformation in ABM (A) Rule map: angle to the gradient (x-axis) vs traction strength (y-axis). DD traction is strongest orthogonal to Wnt5a; DP is strongest along Fgf8b; PP is weak/constant. Dotted boxes show traction biases in DD and DP cases. Bottom box shows Tissue-scale CE: Distal motion orthogonal to the Wnt5a gradient drives elongation. (B) ABM combining differential adhesion (DD>DP>PP), DP bias toward Fgf8b, Distal Wnt5a production, and DD bias orthogonal to Wnt5a. Dotted line marks protrusion shape. (C) CE vectors of ABM simulations projected onto the elongation axis. (D) Schematic for the forces in the ABM: neighbors within diameter R exert identity-/gradient-dependent tractions; short-range repulsion prevents overlap; labels link terms to the equations. (E) Evolution of Distal cell position: sum of DD/DP tractions with the above angle biases plus repulsion within δ. (F) Evolution of Proximal cell position: sum of PP/DP tractions plus repulsion within δ. (G) Evolution of Wnt5a concentration: Distal production, slow diffusion, and decay. (H) Evolution of Fgf8b concentration: slow diffusion with decay inside organoids and a saturated external reservoir, forming an inner–outer gradient. (I) Overview of gradients in G–H. Next, we further validated our model by directly comparing the observed cellular events in organoids with ABM predictions. Consistent with experimental observations, the maximal convergent extension was detected at the distal-proximal cell interface in the ABM ( Fig. 3E , 4C ). In contrast, when only differential adhesion and Fgf8b-directed traction bias are implemented ( Fig. 2E , Supplementary movie 5), the convergent extension maximum is weaker and co-localizes with the distal cluster ( Fig. S7B ). The ABM also recapitulated the increase in the number of protrusions as the initial cell number increased in the organoids ( Fig. S8A , top row). Furthermore, the ABM predicted that as the proximal cell ratio increases, the number of protrusions also increases ( Fig. S8A , bottom row). Experimental data showed that when the distal cell ratio is higher, distal cell clusters do not bifurcate by themselves ( Fig. S8B, C ). In contrast, when proximal cell ratios were higher, it increased the chances of distal cluster separation ( Fig S8D ), indicating that our ABM is capable of predicting accurately the morphogenic events occurring in our limb-mesenchymal organoids. Download figure Open in new tab Figure S8. Parameters regulating protrusion number. (A) ABM parameters controlling protrusions. Predicted protrusion number varies with (i) total cell number ( N ) and (ii) Distal cell fraction ( D/(D+P) ). The panel summarizes how changes in these parameters shift the distribution of protrusions. (B-D) Starting composition vs protrusions. Left: Stacked bars show the percentage of organoids with 0, 1, 2, or 3 protrusions for each starting condition (D_P). Numbers inside segments = raw counts; x-axis labels include total n; y-axis = %. Colors (0–3) are consistent across panels. Right: Representative images at 72 h; each image is annotated with the starting condition and D:P ratio. (B) Group 1 (P fixed at 6,000). (C) Group 2 (P fixed at 12,000). (D) Group 3 (P increased). The final form of our ABM is a dynamical system defined by four constitutive equations describing the evolution of distal and proximal cell positions, and the concentrations of morphogens for convergent extension (Wnt5a, in this case), and for chemotaxis (Fgf8b, in this case). The first equation indicates that the time evolution of a distal cell position can be determined by the sum of two kinds of traction forces and one kind of repulsion forces ( Fig. 4D, E ). The first kind of traction forces is derived from a sum of elementary forces between distal-distal cell pairs within the distance R . Each of these elementary forces is directed along the vector linking the two cells with a magnitude given by the strong adhesion ( ξ DD : based on the doublet experiments) and traction bias in the direction orthogonal to Wnt5a gradient (for convergent extension). Likewise, distal-proximal traction forces can be determined by the sum of the vectors linking each distal cell with the proximal cells within the distance R of it, weighted by the moderate adhesion ( ξ DP ) and with traction bias along Fgf8b gradient (for chemotaxis). The cell-cell repulsive forces come from the requirement of volume exclusion within the distance δ ( R > δ ) ( Fig. 4D, E ). The time evolution of the proximal cell position can be computed similarly, with the difference that the proximal-proximal traction forces are weighted by weak adhesion ( ξ PP ) and with no influence from neither Fgf8b or Wnt5 gradient ( Fig. 4F ). The time evolution of the morphogen concentrations can be defined by a simple reaction-diffusion network. The evolution of Wnt5a concentration is determined by the decay rate and the production rate, which depends on the distal cell density, and the slow diffusion within the organoid ( Fig. 4G ). Whereas, the evolution of Fgf8b concentration is determined by the decay rate, the saturated diffusion outside the organoid, and the slow diffusion within the organoid ( Fig. 4H, I ). Derivation of a PDE model that exhibits the same mathematical structure as the classical fingering instability As our ABM faithfully recapitulated the key events occurring in our limb-mesenchymal organoids at the cellular-level, we next established a more fundamental continuum framework for the emergence of digit morphologies at the tissue-level. Continuum equations are composed of principled mathematical building blocks amenable to analyzing and interpreting the general dynamics. We started from equations governing the microscopic scale of distal and proximal cell positions. We derived a set of partial differential equations (PDEs) for distal and proximal densities via a mean-field limit 35 when the number of cells is taken to infinity ( N →∞), followed by an appropriate space-time rescaling that makes the interaction range ( R ≡ ε ) much smaller than the size of tissues ( V = 1) to observe the consistent development of patterns at a macroscopic scale. ( Fig. 5A–C ; see Supplementary_MathematicalTheory for details). In this limit, the organoid is described solely by two fields, ρ D and ρ P , interpretable as the probabilities (densities) of finding distal and proximal cells at a given location. As the regions of the non-zero densities define the tissue shape, the resulting PDEs describe the evolution of tissue geometry. The force terms and parameters specified in the ABM carry through to the PDE system, while the morphogen concentrations can be described using the same kind of reaction-diffusion network as in the ABM. Download figure Open in new tab Figure 5. Fingering instability in PDE models. (A) Derivation scheme. The ABM is coarse-grained by letting the cell number N→∞ while keeping the interaction range R constant. At the tissue scale, R ≡ ε is very small relative to tissue size V. ρ ε : Distal and Proximal cell densities. f, w : chemical gradient (Fgf8b and Wnt5a). (B) ABM: Distal/Proximal positions evolve by pairwise tractions plus short-range repulsion. (C) Distal ρ ε,D and Proximal ρ ε,P densities evolution with traction (dark/light blue) and repulsion (red) effects. (D) Advection terms dominate at order ε −1 : Distal density drifts up ν, opposite to Proximal density (see cartoon). (E) At the next order, a Cahn–Hilliard (CH) equation describes the evolution of the Distal cell density: the CH term is coupled with a novel anisotropic diffusion operator arising from biased tractions in the direction w . (F) PDE simulation of digit organoids: endpoint of the simulation of model in panel D with traction bias η =1 and differential adhesion α =4 (Left, see Supplementary and Fig. S9A ), coupled with morphogens evolution following the equations in Fig. 4H and 4I (Right). (G) Same as (F) in a in vivo -like setup initiated from a hemispherical Distal tissue and a constant unidirectional upward ν (see Fig. S9C ). Download figure Open in new tab Figure S9. PDE simulations with variable parameters. (A) Simulations of equations from Fig 5D coupled with the model in Fig. 4H and 4I with traction bias η ∈{0, 0.33, 0.66, 1} and differential adhesion α ∈{1, 2, 3, 4} (see Fig 5G ). Setting identical to the ABM version ( Fig. S7B ). Simulations are stopped when instabilities have developed or at the arbitrary endpoint T =42. (B) Simulations of the Distal tissue only as stated in Fig. 5E initiated from a sphere, with constant radial ν and various traction bias values η and tissue sizes N . Simulations are stopped when instabilities have developed or at the arbitrary endpoint T =10. (C) Same as Fig 5G for various traction bias values η and tissue sizes N . (D) Normalized value c k / c 0 of the largest Fourier coefficient of each contour as a function of η in the setting of Fig S9B . Instabilities develop if c k / c 0 reaches the arbitrary cap value 0.05. (E) In the setting of Fig S9B, heatmap of the values of the first 12 Fourier modes of the contours as a function of N for η = 0.6 at the endpoint of the simulation. The number of protrusions corresponds to the maximal Fourier mode, normalized to 1. The resulting system is composed of three fundamental components ( Fig. 5D–E ). First, an advection term which acts on a fast timescale (the coefficient ε −1 , which is large as ε is small). This term corresponds to Fgf8b-driven chemotaxis and produces rapid cell sorting, as distal and proximal cells migrate in opposite directions along the gradient ( Fig. 5D ). The remaining terms —arising from differential adhesion and Wnt5a-directed traction— take the form of a Cahn–Hilliard equation 36 ( Fig. 5E ), which in a similar limit was also observed for other multicellular systems 37 – 39 . Here, we have extended this latter approach to two cell populations coupled with a reaction-diffusion network of morphogens.The Cahn–Hilliard term stabilizes the density towards 0 or 1, corresponding respectively to the outside region (medium) and the tissue. If there were only the Cahn–Hilliard term, the tissue would be stable as a sphere. However, the traction bias introduces a novel anisotropic diffusion term, whose directionality is set by the Wnt5a gradient. This anisotropic diffusion term alone would uniformly expand the tissue and dilute its own density. Here, the combination of the Cahn–Hilliard term (as a stabilizing factor) and the anisotropic diffusion (as an expanding factor) makes the interface unstable and creates protrusions ( Fig. 5F, S9A, B ), a process mathematically linked to the classical "fingering instability" in fluid dynamics (see Discussion) 36 , 40 . As in the ABM, the number of protrusions is controlled by the strength of the traction bias (η), tissue size (N), and the ratio of distal cells ( Fig. S9B–E ). In a setup mimicking limb development in vivo —where morphogenesis initiates from a hemispherical distal tissue with a unidirectional Fgf8b/Wnt5a gradient— digit-like protrusions emerge in the distal direction ( Fig. 5G , S9C). Thus, although the distal tissue is not a viscous fluid per se, its shape evolution is well described by the same class of PDEs, appearing as finger-like protrusions. Despite its name, fingering instability has not been applied to-date as a mechanistic explanation for the formation of digits or “fingers”. Nevertheless, our work shows that fingering instability is sufficient to explain the emergence of digit-like structures. Discussion To establish an explanatory framework for the emergence of digit tissue-like morphologies, we developed a novel limb mesenchymal organoid capable of forming cell aggregates that spontaneously break morphological symmetry to elongate digit-like structures. Gene expression profiling and molecular characterization showed that our organoids recapitulate many of the key cellular and molecular mechanisms of digit formation during development in vivo . Based on our bottom-up approach from an ABM and biological validation, we demonstrate that the tissue-level symmetry-breaking and elongation can be achieved by: 1) an inner-outer gradient of Fgf8b, 2) differential adhesion, 3) a self-organizing Wnt5a gradient, 4) chemotaxis, and 5) convergent-extension. Furthermore, by deriving a continuum description, our PDE derivation reveals that fingering instability is an underlying analytical feature of this ABM, which appears to be the key mechanism driving the emergence of digit morphologies, drawing an unexpected mathematical parallel between tissue morphogenesis and fluid physics. Interestingly, by varying the total cell number and proportion of distal/proximal cells in our model, we were able to mimic phenotypes previously reported in Hox and Gli3 mutants 5 . In Gli3 mutants, limb buds enlarge and digit number increases—consistent with our observation that larger organoids produce more digit-like protrusions ( Fig. 1B, C , S8A, top row). While increasing the proportion of distal cells does not enhance the bifurcation of distal clusters, increasing proximal cells yields a higher number of digits in both organoids and our simulations ( Fig. S8A , bottom row, B-D). In mutant embryos with complete knock-out of distal Hox genes, digits were abolished, similar to the absence of distal cells in our organoids. In contrast, a partial knock-out of distal Hox genes increases the number of digits, which is consistent with our organoid data containing fewer distal cells generating more digits ( Fig. S8 ). Fingering instability is a mechanochemical morphogenesis process that differs fundamentally from the classical Turing instabilities. During normal limb development in vivo , distal cells are typically not in direct contact with proximal cells and symmetry-breaking has —at least at the molecular patterning level— been explained by Turing instabilities driven by biochemical reaction-diffusion across the anterior-posterior axis 5 , 6 . Previous 2D cultures of distal cells appear to capture this early phase of molecular pattern symmetry-breaking 22 , 5 , 6 .( Fig S4B, C ). In contrast, our experimental approach —mixing distal and proximal cells in 2D/3D— allows us to capture morphological symmetry-breaking and subsequent digit elongation along the distal-proximal axis. Although the distal-proximal cell fate of the limb bud is specified by the biochemical gradients of Fgf8 and Wnt3 14 , 15 , cell sorting also likely contributes to the establishment of proximal-distal domains during limb development in vivo , as observed in our model 24 , 41 . Furthermore, the elongation of distal clusters via a Wnt5a-dependent gradient, inducing convergent extension in our organoids, is also well aligned with the mechanisms of digit elongation shown during limb development in vivo 27 – 29 . Therefore, our model not only closely reflects the morphogenetic mechanisms of distal autopod tissue establishment but also shows that digit elongation along the distal-proximal axis occurs based on mechanochemical instability ( Fig. S10 ). Furthermore, the progressive accumulation of distal cells coincided with cartilage differentiation in the 2D micromass culture ( Fig. S3D, E ), although cell-cell biochemical interactions were not extensively analysed. Together, these observations suggest that patterning and morphogenesis are highly coordinated, likely utilizing a combination of both mechanical and biochemical interactions. Depending on the context (e.g. stage and species), one mechanism may dominate over the other. In mouse embryos, the earliest stage of digit patterning involves a periodic arrangement of molecular signals, which is seen hours before any morphological changes are found 6 , whereas amphibians form their morphological digits sequentially, a process where mechanical interactions likely dominate 42 . Download figure Open in new tab Figure S10. Biological interpretation of the existing and current model. Condition A in Fig S4 corresponds to establishment of the distal autopod domain and elongation of digits along the proximal–distal (PD) axis. Conditions B and C in Fig S4 correspond to the periodic patterning of digit and interdigital tissues along the anterior–posterior axis at, respectively, the distal-most part (Condition B: longer wavelength) and a more proximal region (Condition C) of the distal autopod. Fingering instability is often discussed in the context of the Saffman–Taylor problem 36 , 40 , where a viscous fluid invades into a less viscous one. In general, the fluid interface tends to be round or flat due to the surface tension (stability), but when coupled to an expansive drive, the interface may develop protrusions (fingering). The key here is a mathematical balance between stabilization and expansion factors, which are described by the Cahn-Hilliard equation coupled with additional advection terms. While we did not set out to equate digit morphogenesis with the physics of the Saffman-Taylor instability, the same mathematical framework arose unexpectedly from the continuum limit of our agent-based model. We would like to stress that fingering instabilities should be considered as a general emergent phenomenon triggered by the balance between stabilizing and expansive factors and not necessarily linked to viscosity effects. Here, although of different physical nature, anisotropic diffusion plays the same mathematical role as classical pressure and fluid advection mechanisms in triggering fingering instabilities. To the best of our knowledge, fingering instability based on the Cahn-Hilliard equation, has only been applied to explain biological phenomena such as tumour growth and expansion of epithelial colonies 43 , 44 . Hence, our study is the first to describe digit or “finger” formation using a fingering instability mechanism, where cell sorting and traction bias play a distinctive role to control morphological changes. It could be expected that this general mechanism may also appear in other developmental contexts where tissue deformation is driven by convergent extension such as body-axis elongation. Recently, a growing number of studies leveraging multicellular in vitro models have revealed mechanisms underlying emergent tissue morphogenesis. In epithelial tissues, such as the optic cup and the intestinal crypt, morphogenesis can be described as a tissue sheet deformation and the interactions between the mechanical properties of the tissues and the cell fates can explain self-organizing morphogenesis 11 , 45 – 48 . These insights have enabled robust and precise control in tissue engineering 49 , 50 . With mesenchymal tissues, “gastruloids” that model the development of growing axial tissues, show similar morphological characteristics as our model, including symmetry-breaking and elongation, that are based on differential adhesions and convergent extension 51 – 53 . Hence, these studies with ours, provide insights into the previously unappreciated dynamic relationship between the genotypic variation controlling microscopic cellular behaviours and the macroscopic tissue-level morphological phenotypes. Lastly, with our bottom-up approach, since the parameters of our PDE are directly linked to our ABM, the parameters governing the tissue-level morphology can, in principle, be estimated empirically. This may also reveal underlying mechanisms governing evolutional conservation or variation of the digit number. For instance, despite the size of limb buds varying across animal species and number of digits being determined by the size of the limb bud, the five digit phenotype is robustly conserved among most extant tetrapods 54 . While some animals exhibit digit reduction via ad hoc developmental mechanisms 55 , 56 , virtually no extant species has more than five digits. The empirical mechanisms of this strangely conserved phenotype, often referred to as developmental constraint, might be explained by rigorously characterizing the relationship between the microscopic developmental parameters and the macroscopic tissue morphologies. Overall, our work offers a conceptual framework which paves a path for uncovering previously unexplained causes of phenotypic variations in tissue morphologies and for precision tissue engineering strategies. Supplementary movie 1. Limb-mesenchymal organoid: culture conditions. Limb mesenchyme from Sox9–KI–EGFP E11.5 mouse embryos is aggregated and cultured for 193 h in NC or FW (Fgf8+Wnt3a) conditions, initiated with 9,000 or 27,000 cells; Sox9-GFP marks cartilage (green) on bright-field. Associated with Fig 1B . Supplementary movie 2. Limb-mesenchymal organoid: the distal:proximal labelling. Distal–proximal lineage labeling. The cells from the distal half of autopod are labeled H2B-Venus (green) and proximal cells H2B-tdTomato (magenta) before aggregation, and cultured in an organoid for 76 h in NC and FW, started from 9,000 or 27,000 cells, overlaid with bright-field. Tips were enriched in Distal cells. Associated with Fig. 1E . Supplementary movie 3. 2D micromass: conditions. The cells from the distal and proximal half of autopod of Sox9–KI–EGFP (white) E11.5 mouse embryos are labeled before aggregation, and cultured as micromass culture. (B) Micromass culture of limb mesenchyme in NC and FW at 72 h; Sox9-GFP in white; Distal and Proximal cells labeled with H2BtdTomato (green) and H2BtagBFP (magenta). Associated with Fig. S3B . Supplementary movie 4. Simulation: differential adhesions. Agent-based model implementing differential adhesion only (DD>DP>PP). Right: Fgf8b and Wnt5a fields. Associated with Fig. 2B . Supplementary movie 5. Simulation: differential adhesions + Fgf8b-driven chemotaxis. Agent-based model implementing differential adhesion + traction bias toward the Fgf8b gradient. Right: Fgf8b and Wnt5a fields. Associated with Fig. 2E . Supplementary movie 6. Limb-mesenchymal organoid: merging with chemoattractant aggregate. A proximal signaling aggregate (tagBFP or Fgf8b-P2A-tagBFP, blue) is assembled with a Distal(Green)+Proximal(Magenta) aggregate and cultured in NC condition for 96 h. Associated with Fig. 2G and S5D. Supplementary movie 7. Limb-mesenchymal organoid: overexpression of the distal morphogens. Overexpression of tagBFP or Wnt5a-P2A-tagBFP (blue) in Proximal cells, cultured in FW condition for 96h. Associated with Fig. 3F . Supplementary movie 8. Simulation: differential adhesions + Fgf8b-driven chemotaxis + Wnt5a-driven convergent extension. ABM combining differential adhesion (DD>DP>PP), DP bias toward Fgf8b, Distal Wnt5a production, and DD bias orthogonal to Wnt5a. Right: Fgf8b and Wnt5a fields. Associated with Fig. 4B . Funding This work was supported by ASHBi Fusion Research Grant (FY2023-2024: to R.T, S.P., and A.D.), Japan Society for the Promotion of Science (JSPS) Grant-in-Aid for Research Activity Start-up (JP20K22652 to R.T.), Grant-in-Aid for Early-Career Scientists (JP22K15122 to R.T.), KAKENHI Grant-in-Aid for Early-Career Scientists (JP23K13015 to A.D., JP24K16962 to S.P.), AND PRESTO, Japan Science and Technology Agency (JST) (JPMJPR2537 to R.T.). Grant-in-Aid for Transformative Research Areas (Ministry of Education, Culture, Sports, Science, and Technology, MEXT), Japan (JP23H04933 to M.E., JP22H05110 supporting S.P.), and Core Research for Evolutional Science and Technology (CREST, JST) (JPMJCR12W2 to M.E.). Author contributions R.T. conceived, designed the study and performed experiments. A.D. and S.P. performed simulations and mathematical analyses. R.K. and S.O. performed PIC experiments. K.T. performed doublet analysis. R.N. performed scRNA-seq multiplexing. H.A. provided the Sox9-GFP mouse line. R.Takada and S.T. produced the Wnt5a antibody. Cauchy–Green tensor analysis: S.P. and M.M., with initial analyses by M.M., supervised by J.S., who also provided critical scientific feedback and discussions. M.E. conceived, designed and supervised the study. R.T., A.D., and S.P. wrote the manuscript with input from all authors. Competing interests Authors declare no conflict of interests. Materials & Correspondence R.T. ( tsutsumi.rio.8x{at}kyoto-u.ac.jp ) and M.E. ( eiraku.mototsugu.2u{at}kyoto-u.ac.jp ). Supplementary Contents Supplementary_scRNA-seq Supplementary_MathematicalTheory Supplementary_CauchyGreenTensorEstimation Supplementary movie 1 (Limb-mesenchymal organoid: culture conditions) Supplementary movie 2 (Limb-mesenchymal organoid: the distal:proximal labelling) Supplementary movie 3 (2D micromass: conditions) Supplementary movie 4 (Simulation: differential addhesions) Supplementary movie 5 (Simulation: differential adhesions+Fgf8b-driven chemotaxis) Supplementary movie 6 (Limb-mesenchymal organoid: merging with chemoattractant aggregate) Supplementary movie 7 (Limb-mesenchymal organoid: overexpression of the distal morphogens) Supplementary movie 8 (Simulation: differential adhesions+Fgf8b-driven chemotaxis+Wnt5a-driven convergent extension) Acknowledgements We thank all members of the Eiraku Lab for helpful discussions; Dr. Hitomi Watanabe for generating Sox9–GFP mouse embryos; and Dr. Ko Sugawara for technical support with the ELEPHANT plug-in; and Dr. Elazer Zelzer for sharing scRNAseq data. We thank Single-cell Genome Information Analysis Core (SignAC) at WPI-ASHBi (Kyoto University) for their support. We also thank Dr. Spyros Goulas at WPI-ASHBi (Kyoto University) for critical reading and constructive suggestions on our manuscript. Funder Information Declared Japan Society for the Promotion of Science , JP20K22652 , JP22K15122 , JP23K13015 , JP24K16962 , JP23H04933 Japan Science and Technology Agency , JPMJPR2537 , JPMJCR12W2 ASHBi , Fusion grant Footnotes 1. The clarification on the relationship between our framework and Turing models that has been previously published. 2. Addition of the experiments showing the distal cells exihibit chemotaxis towards Fgf8b. 3. Edits for improving visibility. References 1. ↵ Saxena , A. , Towers , M. & Cooper , K. L . The origins, scaling and loss of tetrapod digits . Phil Trans R Soc B 372 , 20150482 ( 2017 ). OpenUrl CrossRef PubMed 2. ↵ Martinez Arias , A. The Master Builder: How the New Science of the Cell Is Rewriting the Story of Life . ( Basic Books , New York, NY , 2023 ). 3. ↵ Anderson , P. W. More Is Different . Science 177 , 393 – 396 ( 1972 ). OpenUrl FREE Full Text 4. ↵ Turing , A. M . The chemical basis of morphogenesis . Philos. Trans. R. Soc. Lond. B Biol. Sci . 237 , 37 – 72 ( 1952 ). OpenUrl CrossRef 5. ↵ Sheth , R. et al. Hox Genes Regulate Digit Patterning by Controlling the Wavelength of a Turing-Type Mechanism . Science 338 , 1476 – 1480 ( 2012 ). OpenUrl Abstract / FREE Full Text 6. ↵ Raspopovic , J. , Marcon , L. , Russo , L. & Sharpe , J . Digit patterning is controlled by a Bmp-Sox9-Wnt Turing network modulated by morphogen gradients . Science 345 , 566 – 570 ( 2014 ). OpenUrl Abstract / FREE Full Text 7. ↵ Green , J. B. A . Resolving morphogenesis into quantifiable cell behaviours . Development 149 , dev199794 ( 2022 ). 8. ↵ Sasai , Y . Cytosystems dynamics in self-organization of tissue architecture . Nature 493 , 318 – 326 ( 2013 ). OpenUrl CrossRef PubMed Web of Science 9. Warmflash , A. , Sorre , B. , Etoc , F. , Siggia , E. D. & Brivanlou , A. H. A method to recapitulate early embryonic spatial patterning in human embryonic stem cells . Nat. Methods 11 , 847 – 854 ( 2014 ). OpenUrl CrossRef PubMed 10. Eiraku , M. et al. Self-Organized Formation of Polarized Cortical Tissues from ESCs and Its Active Manipulation by Extrinsic Signals . Cell Stem Cell 3 , 519 – 532 ( 2008 ). OpenUrl CrossRef PubMed Web of Science 11. ↵ Eiraku , M. et al. Self-organizing optic-cup morphogenesis in three-dimensional culture . Nature 472 , 51 ( 2011 ). OpenUrl CrossRef PubMed Web of Science 12. ↵ Kumar , D. & Lassar , A. B . Fibroblast Growth Factor Maintains Chondrogenic Potential of Limb Bud Mesenchymal Cells by Modulating DNMT3A Recruitment . Cell Rep . 8 , 1419 – 1431 ( 2014 ). OpenUrl CrossRef PubMed 13. ↵ ten Berge , D. , Brugmann , S. A. , Helms , J. A. & Nusse , R . Wnt and FGF signals interact to coordinate growth with cell fate specification during limb development . Development 135 , 3247 – 3257 ( 2008 ). OpenUrl Abstract / FREE Full Text 14. ↵ Cooper , K. L. et al. Initiation of Proximal-Distal Patterning in the Vertebrate Limb by Signals and Growth . Science 332 , 1083 – 1086 ( 2011 ). OpenUrl Abstract / FREE Full Text 15. ↵ Roselló-Díez , A. , Ros , M. A. & Torres , M. Diffusible Signals, Not Autonomous Mechanisms, Determine the Main Proximodistal Limb Subdivision . Science 332 , 1086 – 1088 ( 2011 ). OpenUrl Abstract / FREE Full Text 16. ↵ Honda , M. et al. High-depth spatial transcriptome analysis by photo-isolation chemistry . Nat. Commun . 12 , 4416 ( 2021 ). OpenUrl CrossRef PubMed 17. ↵ Bruneau , S. , Johnson , K. R. , Yamamoto , M. , Kuroiwa , A. & Duboule , D . The Mouse Hoxd13 spdh Mutation, a Polyalanine Expansion Similar to Human Type II Synpolydactyly (SPD), Disrupts the Function but Not the Expression of Other Hoxd Genes . Dev. Biol . 237 , 345 – 353 ( 2001 ). OpenUrl CrossRef PubMed Web of Science 18. Houzelstein , D. , Cohen , A. , Buckingham , M. E. & Robert , B . Insertional mutation of the mouse Msx1 homeobox gene by an nlacZ reporter gene . Mech. Dev . 65 , 123 – 133 ( 1997 ). OpenUrl CrossRef PubMed Web of Science 19. Tzchori , I. et al. LIM homeobox transcription factors integrate signaling events that control three-dimensional limb patterning and growth . Development 136 , 1375 – 1385 ( 2009 ). OpenUrl Abstract / FREE Full Text 20. ↵ Hargrave , M. et al. Expression of the Sox11 gene in mouse embryos suggests roles in neuronal maturation and epithelio-mesenchymal induction . Dev. Dyn . 210 , 79 – 86 ( 1997 ). OpenUrl CrossRef PubMed Web of Science 21. ↵ Summerbell , D. , Lewis , J. H. & Wolpert , L . Positional Information in Chick Limb Morphogenesis . Nature 244 , 492 – 496 ( 1973 ). OpenUrl CrossRef PubMed Web of Science 22. ↵ Fuiten , A. M. , Yoshimoto , Y. , Shukunami , C. & Stadler , H. S . Digits in a dish: An in vitro system to assess the molecular genetics of hand/foot development at single-cell resolution . Front. Cell Dev. Biol . 11 , ( 2023 ). 23. ↵ Wada , N . Spatiotemporal changes in cell adhesiveness during vertebrate limb morphogenesis . Dev. Dyn . 240 , 969 – 978 ( 2011 ). OpenUrl CrossRef PubMed 24. ↵ Stadler , H. S. , Higgins , K. M. & Capecchi , M. R . Loss of Eph-receptor expression correlates with loss of cell adhesion and chondrogenic capacity in Hoxa13 mutant limbs . Development 128 , 4177 – 4188 ( 2001 ). OpenUrl Abstract / FREE Full Text 25. ↵ Maître , J.-L. et al. Adhesion Functions in Cell Sorting by Mechanically Coupling the Cortices of Adhering Cells . Science 338 , 253 – 256 ( 2012 ). OpenUrl Abstract / FREE Full Text 26. ↵ Lesnicar-Pucko , G. , Belmonte , J. M. , Musy , M. , Glazier , J. A. & Sharpe , J . Cellular mechanisms of chick limb bud morphogenesis . bioRxiv 2020.09.10.292359 ( 2020 ) doi: 10.1101/2020.09.10.292359 . OpenUrl Abstract / FREE Full Text 27. ↵ Gao , B. et al. Wnt Signaling Gradients Establish Planar Cell Polarity by Inducing Vangl2 Phosphorylation through Ror2 . Dev. Cell 20 , 163 – 176 ( 2011 ). OpenUrl CrossRef PubMed Web of Science 28. Gao , B. et al. Coordinated directional outgrowth and pattern formation by integration of Wnt5a and Fgf signaling in planar cell polarity . Development 145 , dev163824 ( 2018 ). 29. ↵ Parada , C. et al. Mechanical feedback defines organizing centers to drive digit emergence . Dev. Cell 57 , 854 – 866 .e6 ( 2022 ). OpenUrl CrossRef PubMed 30. ↵ Tinevez , J.-Y. et al. TrackMate: An open and extensible platform for single-particle tracking . Methods 115 , 80 – 90 ( 2017 ). OpenUrl CrossRef PubMed 31. ↵ Sugawara , K. , Çevrim , Ç. & Averof , M . Tracking cell lineages in 3D by incremental deep learning . eLife 11 , e69380 ( 2022 ). OpenUrl CrossRef PubMed 32. ↵ Sheth , R. et al. Distal Limb Patterning Requires Modulation of cis -Regulatory Activities by HOX13 . Cell Rep . 17 , 2913 – 2926 ( 2016 ). OpenUrl CrossRef PubMed 33. ↵ Salsi , V. , Vigano , M. A. , Cocchiarella , F. , Mantovani , R. & Zappavigna , V . Hoxd13 binds in vivo and regulates the expression of genes acting in key pathways for early limb and skeletal patterning . Dev. Biol . 317 , 497 – 507 ( 2008 ). OpenUrl CrossRef PubMed Web of Science 34. ↵ Montero , J. A. , Lorda-Diez , C. I. , Gañan , Y. , Macias , D. & Hurle , J. M . Activin/TGFβ and BMP crosstalk determines digit chondrogenesis . Dev. Biol . 321 , 343 – 356 ( 2008 ). OpenUrl CrossRef PubMed Web of Science 35. ↵ Chaintron , L.-P. & Diez , A . Propagation of chaos: A review of models, methods and applications . I. Models and methods. Kinet. Relat. Models 15 , 895 – 1015 ( 2022 ). OpenUrl 36. ↵ Cahn , J. W. & Hilliard , J. E . Free Energy of a Nonuniform System. I. Interfacial Free Energy . J. Chem. Phys . 28 , 258 – 267 ( 1958 ). OpenUrl CrossRef Web of Science 37. ↵ Armstrong , N. J. , Painter , K. J. & Sherratt , J. A . A continuum approach to modelling cell–cell adhesion . J. Theor. Biol . 243 , 98 – 113 ( 2006 ). OpenUrl CrossRef PubMed Web of Science 38. Falcó , C. , Baker , R. E. & Carrillo , J. A . A Local Continuum Model of Cell-Cell Adhesion . SIAM J. Appl. Math . 84 , S17 – S42 ( 2024 ). OpenUrl 39. ↵ Kim , I. , Mellet , A. & Sheung-Him Wu , J . Aggregation-Diffusion Phenomena: From Microscopic Models to Free Boundary Problems . in Active Particles , Volume 4 217 – 281 ( Birkhäuser, Cham , 2024 ). doi: 10.1007/978-3-031-73423-6_6 . OpenUrl CrossRef 40. ↵ Saffman , P. G. & Taylor , G. I. The penetration of a fluid into a porous medium or Hele-Shaw cell containing a more viscous liquid . Proc. R. Soc. Lond. Ser. Math. Phys. Sci . 245 , 312 – 329 ( 1997 ). OpenUrl 41. ↵ Mizuno , K. , Hirashima , T. & Toda , S . Robust tissue pattern formation by coupling morphogen signal and cell adhesion . EMBO Rep . 25 , 4803 – 4826 ( 2024 ). OpenUrl PubMed 42. ↵ Fröbisch , N. B. & Shubin , N. H . Salamander limb development: Integrating genes, morphology, and fossils . Dev. Dyn . 240 , 1087 – 1099 ( 2011 ). OpenUrl CrossRef PubMed 43. ↵ Ebenbeck , M. , Garcke , H. & Nürnberg , R . Cahn–Hilliard–Brinkman systems for tumour growth . Discrete Contin. Dyn. Syst. - S 14 , 3989 – 4033 ( 2021 ). OpenUrl 44. ↵ Ye , Y. & Lin , J . Fingering Instability Accelerates Population Growth of a Proliferating Cell Collective . Phys. Rev. Lett . 132 , 018402 ( 2024 ). OpenUrl PubMed 45. ↵ Okuda , S. et al. Strain-triggered mechanical feedback in self-organizing optic-cup morphogenesis . Sci. Adv . 4 , eaau1354 ( 2018 ). 46. Pérez-González , C. et al. Mechanical compartmentalization of the intestinal organoid enables crypt folding and collective cell migration . Nat. Cell Biol . 23 , 745 – 757 ( 2021 ). OpenUrl CrossRef PubMed 47. Yang , Q. et al. Cell fate coordinates mechano-osmotic forces in intestinal crypt formation . Nat. Cell Biol . 23 , 733 – 744 ( 2021 ). OpenUrl CrossRef PubMed 48. ↵ Serra , D. et al. Self-organization and symmetry breaking in intestinal organoid development . Nature 569 , 66 – 72 ( 2019 ). OpenUrl CrossRef PubMed 49. ↵ Gjorevski , N. et al. Tissue geometry drives deterministic organoid patterning . Science 375 , eaaw9021 ( 2022 ). 50. ↵ Nikolaev , M. et al. Homeostatic mini-intestines through scaffold-guided organoid morphogenesis . Nature 585 , 574 – 578 ( 2020 ). OpenUrl CrossRef PubMed 51. ↵ Brink , S. C. van den et al. Symmetry breaking, germ layer specification and axial organisation in aggregates of mouse embryonic stem cells . Development 141 , 4231 – 4242 ( 2014 ). OpenUrl Abstract / FREE Full Text 52. Jong , M. A. de et al. The shapes of elongating gastruloids are consistent with convergent extension driven by a combination of active cell crawling and differential adhesion . PLOS Comput. Biol . 20 , e1011825 ( 2024 ). OpenUrl CrossRef PubMed 53. ↵ McNamara , H. M. , Solley , S. C. , Adamson , B. , Chan , M. M. & Toettcher , J. E . Recording morphogen signals reveals mechanisms underlying gastruloid symmetry breaking . Nat. Cell Biol . 26 , 1832 – 1844 ( 2024 ). OpenUrl CrossRef PubMed 54. ↵ Galis , F. , J.M. van Alphen , J. & A. J. Metz , J. Why five fingers? Evolutionary constraints on digit numbers . Trends Ecol. Evol . 16 , 637 – 646 ( 2001 ). OpenUrl CrossRef Web of Science 55. ↵ Cooper , K. L. et al. Patterning and post-patterning modes of evolutionary digit loss in mammals . Nature 511 , 41 – 45 ( 2014 ). OpenUrl CrossRef GeoRef PubMed 56. ↵ Lopez-Rios , J. et al. Attenuated sensing of SHH by Ptch1 underlies evolution of bovine limbs . Nature 511 , 46 – 51 ( 2014 ). OpenUrl CrossRef PubMed 57. Lewandowski , J. P. , Pursell , T. A. , Rabinowitz , A. H. & Vokes , S. A . Manipulating gene expression and signaling activity in cultured mouse limb bud cells . Dev. Dyn . 243 , 928 – 936 ( 2014 ). OpenUrl CrossRef PubMed 58. Ershov , D. et al. TrackMate 7: integrating state-of-the-art segmentation algorithms into tracking pipelines . Nat. Methods 19 , 829 – 832 ( 2022 ). OpenUrl CrossRef PubMed 59. Tsutsumi, R. RioTsutsumi0617/Yokogawa-CV1000-image-processing-: This is to process multi-well, tile, time, and channel tiff images taken by Yokogawa-CV1000 into simpler view using imageJ . https://github.com/RioTsutsumi0617/Yokogawa-CV1000-image-processing-. 60. mastodon-sc/mastodon. Mastodon Science ( 2025 ). 61. Imoto , Y . Comprehensive noise reduction in single-cell data with the RECODE platform. Cell Rep . Methods 5 , ( 2025 ). 62. Bergen , V. , Lange , M. , Peidli , S. , Wolf , F. A. & Theis , F. J . Generalizing RNA velocity to transient cell states through dynamical modeling . Nat. Biotechnol . 38 , 1408 – 1414 ( 2020 ). OpenUrl CrossRef PubMed 63. Street , K. et al. Slingshot: cell lineage and pseudotime inference for single-cell transcriptomics . BMC Genomics 19 , 477 ( 2018 ). OpenUrl CrossRef PubMed 64. ↵ Markman , S. et al. A single-cell census of mouse limb development identifies complex spatiotemporal dynamics of skeleton formation . Dev. Cell ( 2023 ) doi: 10.1016/j.devcel.2023.02.013 . OpenUrl CrossRef PubMed 65. Ruiz-Villalba , A. et al. Reference genes for gene expression studies in the mouse heart . Sci. Rep . 7 , 24 ( 2017 ). OpenUrl CrossRef PubMed 66. Zhang , K. et al. Silencing of Vangl2 attenuates the inflammation promoted by Wnt5a via MAPK and NF-κB pathway in chondrocytes . J. Orthop. Surg. 16 , 136 ( 2021 ). OpenUrl 67. Paszke , A. et al. PyTorch: An Imperative Style, High-Performance Deep Learning Library . in Advances in Neural Information Processing Systems vol. 32 ( Curran Associates, Inc. , 2019 ). 68. Charlier , B. , Feydy , J. , Glaunès , J. A. , Collin , F.-D. & Durif , G . Kernel Operations on the GPU, with Autodiff, without Memory Overflows . J. Mach. Learn. Res . 22 , 1 – 6 ( 2021 ). OpenUrl 69. Diez , A . SiSyPHÈ: A Python package for the Simulation of Systems of interacting mean-field Particles with High Efficiency . J. Open Source Softw . 6 , 3653 ( 2021 ). OpenUrl 70. Baratta , I. A. , et al. DOLFINx: The next generation FEniCS problem solving environment . Preprint at doi: 10.5281/zenodo.10447666 ( 2023 ). OpenUrl CrossRef 71. ParaView: An End-User Tool for Large-Data Visualization . in Visualization Handbook 717 – 731 ( Butterworth-Heinemann , 2005 ). doi: 10.1016/B978-012387582-2/50038-1 . OpenUrl CrossRef View the discussion thread. Back to top Previous Next Posted April 06, 2026. Download PDF Supplementary Material Email Thank you for your interest in spreading the word about bioRxiv. NOTE: Your email address is requested solely to identify you as the sender of this article. Your Email * Your Name * Send To * Enter multiple addresses on separate lines or separate them with commas. You are going to email the following Mechanochemical instabilities drive digit morphogenesis in organoids Message Subject (Your Name) has forwarded a page to you from bioRxiv Message Body (Your Name) thought you would like to see this page from the bioRxiv website. Your Personal Message CAPTCHA This question is for testing whether or not you are a human visitor and to prevent automated spam submissions. Share Mechanochemical instabilities drive digit morphogenesis in organoids Rio Tsutsumi , Antoine Diez , Steffen Plunder , Ryuichi Kimura , Oki Shinya , Kaori Takizawa , Rei Nakano , Haruhiko Akiyama , Ritsuko Takada , Shinji Takada , Marco Musy , James Sharpe , Mototsugu Eiraku bioRxiv 2025.08.31.673315; doi: https://doi.org/10.1101/2025.08.31.673315 Share This Article: Copy Citation Tools Mechanochemical instabilities drive digit morphogenesis in organoids Rio Tsutsumi , Antoine Diez , Steffen Plunder , Ryuichi Kimura , Oki Shinya , Kaori Takizawa , Rei Nakano , Haruhiko Akiyama , Ritsuko Takada , Shinji Takada , Marco Musy , James Sharpe , Mototsugu Eiraku bioRxiv 2025.08.31.673315; doi: https://doi.org/10.1101/2025.08.31.673315 Citation Manager Formats BibTeX Bookends EasyBib EndNote (tagged) EndNote 8 (xml) Medlars Mendeley Papers RefWorks Tagged Ref Manager RIS Zotero Tweet Widget Facebook Like Google Plus One Subject Area Developmental Biology Subject Areas All Articles Animal Behavior and Cognition (7635) Biochemistry (17697) Bioengineering (13895) Bioinformatics (41951) Biophysics (21456) Cancer Biology (18594) Cell Biology (25520) Clinical Trials (138) Developmental Biology (13381) Ecology (19903) Epidemiology (2067) Evolutionary Biology (24323) Genetics (15612) Genomics (22510) Immunology (17738) Microbiology (40401) Molecular Biology (17184) Neuroscience (88622) Paleontology (667) Pathology (2833) Pharmacology and Toxicology (4825) Physiology (7644) Plant Biology (15158) Scientific Communication and Education (2046) Synthetic Biology (4296) Systems Biology (9825) Zoology (2271)