# Pseudodynamic analysis of heart tube formation in the mouse reveals strong regional variability and early left-right asymmetry

Isaac Esteban, Patrick Schmidt, Audrey Desgrange, Morena Raiola, Susana Temiño, Sigolène Meilhac, Leif Kobbelt, Miguel Torres

<figure><img src="image_1.png" /></figure>

## To cite this version:

Isaac Esteban, Patrick Schmidt, Audrey Desgrange, Morena Raiola, Susana Temiño, et al., Pseudodynamic analysis of heart tube formation in the mouse reveals strong regional variability and early left-right asymmetry. Nature Cardiovascular Research, 2022, 1, pp.504-517. 10.1038/s44161-022-00065-1. pasteur-03851032

HAL Id: pasteur-03851032

https://pasteur.hal.science/pasteur-03851032

Submitted on 14 Nov 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

# Pseudodynamic analysis of heart tube formation in the mouse reveals strong regional variability and early left-right asymmetry

Isaac Esteban¹², Patrick Schmidt D³, Audrey Desgrange⁴, Morena Raiola<fnref n="1" />,², Susana Temiño<fnref n="1" />, Sigolène M. Meilhac⁴, Leif Kobbelt³ and Miguel Torres D<fnref n="1" />✉

Understanding organ morphogenesis requires a precise geometrical description of the tissues involved in the process. The high morphological variability in mammalian embryos hinders the quantitative analysis of organogenesis. In particular, the study of early heart development in mammals remains a challenging problem due to imaging limitations and complexity. Here, we provide a complete morphological description of mammalian heart tube formation based on detailed imaging of a temporally dense collection of mouse embryonic hearts. We develop strategies for morphometric staging and quantification of local morphological variations between specimens. We identify hot spots of regionalized variability and identify *Nodal*-controlled left-right asymmetry of the inflow tracts as the earliest signs of organ left-right asymmetry in the mammalian embryo. Finally, we generate a three-dimensional+t digital model that allows co-representation of data from different sources and provides a framework for the computer modeling of heart tube formation

The mammalian embryo is highly regulative, which implies it can consistently develop stereotyped organs despite initial morphological heterogeneity. Heart tube formation involves complex morphogenetic movements and differentiation patterns, through which form and function are acquired simultaneously<fnref n="1" />. Understanding heart tube morphogenesis would greatly benefit from having access to a spatiotemporal representation of the morphology of the tissues involved in the process. Previous attempts to understand the three-dimensional+time (3D+t) complexity of primitive heart tube formation in amniotes have been based on the segmentation and 3D reconstruction of specimens from tissue sections²⁻⁵. This approach has allowed for the reconstruction of a few discrete stages capturing this process but is limited by a demanding methodology and the distortions of section alignment. Additional efforts have included high-resolution episcopic microscopy and live analysis of the developing chicken or mouse heart; however, this was applied to stages beyond primitive heart tube formation and focused on specific processes of cardiogenesis⁶⁻⁹ or was limited by incomplete imaging of the area of interest and low image resolution¹⁰,¹¹.

Cardiac development starts with the deployment of cardiac mesoderm early in gastrulation¹². Soon after colonizing the antero-lateral region of the embryo, the cardiac mesoderm at the anterior lateral plate undergoes a mesenchymal-to-epithelial transition to produce two layers that enclose the celomic cavity. This results in two independent epithelial sheets: the splanchnic mesoderm, which shares basement membrane with the endoderm, and the somatic mesoderm, which shares basement membrane with the ectoderm. The two mesodermal layers meet medially with the head paraxial mesoderm and laterally at the embryonic-extraembryonic border, forming a closed and continuous 3D surface encasing the anterior celomic cavity (future pericardial cavity).

During embryonic days 7.5-8.5 (E7.5-E8.5) in the mouse, the splanchnic and somatic mesoderm layers progressively move apart from each other, forming the pericardial cavity. This occurs in parallel to embryonic folding and foregut pocket formation, which brings the foregut endoderm inside the embryo. This causes the pericardial cavity to rotate inwards around the left-right axis, driven by the rostral folding of the embryo, while the endoderm invaginates toward the anterior. Following foregut pocket formation, the pericardial cavity is limited at its cranial side by the ectodermal plate and the head mesoderm. Caudally, the cavity limits with the foregut diverticulum in its medial aspect, while the posterior-lateral limits of the cavity extend around the foregut invagination, ending at the somatic and splanchnic mesoderm junction.

The splanchnic mesoderm of the forming pericardial cavity contains the cardiac mesoderm that gives rise to the heart tube. Within the cardiac mesoderm, there are two distinct regions: the first heart field (FHF), which differentiates rapidly to the cardiomyocyte fate, and the second heart field (SHF), which remains highly proliferative and undifferentiated until its precursors are progressively incorporated to the heart at later stages. In the early embryo, the FHF extends anterior laterally, closer to the embryonic-extraembryonic boundary, while the SHF occupies the posterior-medial regions adjacent to the FHF. The initially relatively flat cardiac crescent (CC) starts very active morphogenesis around E7.5, growing and acquiring 3D complexity by folding on itself, detaching from the endoderm and coalescing medially to form the primitive tube, with the venous-to-arterial pole axis in a caudal-to-cranial orientation.

In this work, we aimed to obtain a complete description of tissue anatomy during primitive heart tube morphogenesis in the mouse embryo. To this end, we acquired a collection of high-resolution confocal images from whole specimens at high temporal density and used the images to generate a continuous 3D+t description

<footnote n="1">¹Cardiovascular Regeneration Program, Centro Nacional de Investigaciones Cardiovasculares (CNIC), Madrid, Spain. ²Universidad Politécnica de Madrid, Madrid, Spain. ³Visual Computing Institute, RWTH Aachen University, Aachen, Germany. ⁴Unit of Heart Morphogenesis, Université de Paris, Imagine - Institut Pasteur, INSERM UMR1163, Paris, France. ✉e-mail: mtorres@cnic.es</footnote>
of heart tube morphogenesis. Given that heart tube morphogenesis cannot be uncoupled from the morphological evolution of the surrounding tissues, our model incorporates all pericardial cavity mesoderm, endoderm and the endothelium. Besides obtaining a detailed 3D+t atlas of heart tube formation and its associated tissues, we applied statistical methods to study regional variability and left-right patterning during heart tube formation. The model and knowledge generated provide a strong basis for understanding heart morphogenesis, its variability and its congenital alterations.

## Results

**High-resolution spatiotemporal study of the developing heart.** To generate a high-fidelity 3D+t atlas of early mouse heart development, we imaged 52 fixed and cleared¹³ embryos in cardiac developmental stages ranging from the early CC until heart looping. The collection represents nominal ages from approximately E7.75 to E8.5 (18 h). The estimated approximate temporal density is therefore one specimen every 20 min. The imaged embryos carried *Mesp1*<sup>Cre</sup> and *R26R<sup>MitG</sup>* and *R26R<sup>Tomato</sup>* alleles so that all mesoderm in the head and cardiogenic area was labeled with membrane green fluorescent protein (mGFP) and cytoplasmic Tomato, while the rest of the tissues were labeled only with membrane Tomato. We imaged volumes of ~1 mm in the *x* and *y* dimensions and down to a depth of 677 µm, so that the whole spatial domain of the cardiogenic region and associated tissues was acquired. We used an *x/y* pixel size varying from 0.38 µm to 0.49 µm and a *z* step varying from 0.49 µm to 2.0 µm, depending on the specimen. The varying ranges in resolution were adapted for optimization of acquisition time and minimization of photobleaching with the larger specimens. This strategy provided highly detailed orthogonal optical sections at resolutions beyond the need for 3D segmentation at tissue level (Fig. 1c).

We then semiautomatically segmented tissues of interest in each two-dimensional (2D) confocal slice (Methods). The segmented tissues included the mesoderm of the anterior celomic cavity, discriminating both layers of the lateral plate mesoderm, splanchnic and somatic, up to their fusion with the head paraxial mesoderm (Fig. 1d). Within the splanchnic domain, we used morphological and topological criteria to discriminate the differentiated myocardium from the rest of the splanchnic mesoderm. The cells that were detached from the endoderm and encased, or overlayed, endocardial cells were classified as differentiating cardiomyocytes (Fig. 1d). In addition, we segmented the foregut endoderm, the endocardial lumen and the forming aortae (Fig. 1d and Methods). The 2D segmentations of the different tissues were then used to reconstruct primary 3D images (Fig. 1e), which were subsequently smoothened and transformed into meshes (Fig. 1f and Methods). The meshes constitute closed surfaces describing the interface of each segmented tissue with the surrounding tissues or cavities. In this way, we generated

3D models for the 52 specimens of the collection (an example is shown in Fig. 1g). These 3D models are available for download (see Data availability) and can be visualized using, for example, the open-source application ParaView (www.paraview.org).

Preliminary study of five selected specimens in obvious temporal sequence (Fig. 2a,b) showed how the splanchnic mesoderm folds to form the CC. The early CC extends bilaterally and subdivides the undifferentiated splanchnic mesoderm into medial and distal domains. The medial domain contains the SHF precursors, while the distal domain (recently named juxta-cardiac field) contains precursors of the epicardium and the septum transversum and bears myocardial differentiation capacity¹⁴. As the CC undergoes morphogenesis to generate the cardiac tube, both the juxta-cardiac field and its border with the forming heart tube remain morphologically rather stable (Fig. 2a,b), whereas the medial splanchnic mesoderm (containing the SHF) and its border with the forming heart tube undergo a drastic deformation that is essential for the generation of the outflow tract (OFT), allowing the dorsal closure, which transforms the CC into the primitive heart tube (Fig. 2a,b). This transformation is concomitant with the development of the endocardial cavity, the aortae and their connection through the first branchial arch arteries (Fig. 2c,d).

**Morphometric staging of heart tube formation.** We next aimed the temporal ordering of the specimen collection. The examination of the collection of datasets reflected high morphological heterogeneity of the developing heart within groups of embryos at an apparently similar developmental stage, as judged by somite number, head folds and foregut pocket shapes. The morphological heterogeneity affects both general conformation features, such as proportions, degree of dorsal closure or degree of looping, and local morphological aspects, such as non-stereotyped presence of bulges of different sizes and distribution. To tackle the problem of staging despite this variability, we implemented different morphometric strategies. We first obtained midsurfaces¹⁵ for each of the relevant tissues, consisting of a zero-thickness surface resulting from flattening the original surfaces (Fig. 3a), and defined various landmark distances and curves on such surfaces (Fig. 3b and Extended Data Fig. 1). We next studied the proportions between the lengths of pairs of the reference curves/distances (Fig. 3b) and focused on proportions expected to vary continuously during the studied period (See example in Fig. 3a).

Given that the cardiogenic area evolves from a crescent extending from left to right to a tube extending from cranial to caudal, a general trend during the formation of the heart tube is the increase of the height (cranio-caudal size) to width (left-right size) ratio of the tissues involved. We therefore calculated the proportions between different parameters related to the height-to-width ratio variation (Fig. 3a,b, parameters 1 to 6). We then applied the

**Fig. 1 | Segmentation of high-resolution confocal images reveals precise and detailed morphological representations for understanding heart development.** a, Schematic representation of the workflow for image acquisition. The icon below shows the axes of the system of reference; o/n, overnight; PFA, paraformaldehyde; PBT, phosphate buffer saline with 0.5% triton X-100. b, Confocal raw image reconstruction showing the dimensions of the acquired volumes. Mesodermal tissues show the cytoplasm in red and cell membranes in green; the rest of the tissues show membranes in red. Cell nuclei appear in cyan. Red and blue arrowheads indicate the optical sections shown in c, e, Frontal (left) and sagittal (right) optical sections. Dotted-line boxes highlight the cardiac region. d, Top, frontal optical section at the level of the cardiac inflow tract (IFT) region. The orange line shows the segmented mesoderm of the anterior celomic cavity, including the differentiating myocardium, and the foregut endoderm is represented in blue. The foregut pocket lumen is represented using a blue-dot grid. The endocardial lumen is shown in pink, and the aortic lumen is shown in green. Bottom left, lateral sagittal optical section at the level of one of the aortae. Bottom right, the same optical section showing the subtypes of mesoderm of the pericardial cavity. Myocardium is in red, splanchnic mesoderm is in blue, somatic mesoderm is in green, and head paraxial mesoderm is in purple. e, Dorsal (top) and ventral (bottom) 3D representation of the segmented tissues; mes. mesoderm. f, Pipeline for transforming the segmentation images into smooth, closed meshes. g, Dorsal (top) and ventral (bottom) views of the resulting meshes. The black dotted arrow represents the cut trajectory to achieve genus O topology in the splanchnic shapes (Methods). Abbreviations: a, anterior; p, posterior; r, right; l, left; v, ventral; d, dorsal; HF, head folds; HT, heart tube; ECT, ectoderm; PC, pericardial cavity; RSH, right sinus horn; LSH, left sinus horn; END, endoderm; FP, foregut pocket; AOR, aorta; SECT, surface ectoderm; NECT, neuroectoderm; EL, endocardial lumen; MES, mesoderm; PCM, pericardial cavity mesoderm; FPL, foregut pocket lumen; FGE, foregut endoderm; AOL, aortic lumen; MYO, myocardium; SPL, splanchnic mesoderm; SOM, somatic mesoderm; HPM, head paraxial mesoderm.
k-means¹⁶ clustering algorithm for the determination of groups based on differences in the defined parameters and calculated the average silhouette coefficient (ś; Extended Data Fig. 2), which measures the quality of classifications¹⁷. As we increased the number of groups (k) in the k-means classification, the best ś was obtained for the *d*1/*d*2 proportion (Fig. 3c), *d*1 being the length of the border between the myocardium and SHF and *d*2 being the length of the border between the myocardium and juxta-cardiac field (Fig. 3a). The maximum ś was obtained for a k-means classification into 10 groups. We therefore concluded that *d*1/*d*2 was the best single parameter for classifying the collection and established these

10 groups as a reference for study (Fig. 3f). While the establishment of groups is useful to understand the morphological evolution of the heart, the *d*1/*d*2 parameter is a continuous variable and allows the allocation of single specimens along the heart tube developmental trajectory, independently of grouping. Study of the correlation between *d*1/*d*2 and the rest of the parameters provided a view of how good the different classification methods are with respect to *d*1/*d*2 (Extended Data Fig. 1). The best correlate to *d*1/*d*2 was the height-to-width ratio, which showed a high linear correlation (Extended Data Fig. 3a). By contrast, other parameters either showed a poor correlation or did not show a linear correlation, which resulted in

<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /></figure>
<figure><img src="image_8.png" /><figcaption>**Fig. 2** | Observation of qualitative morphometric aspects of heart tube morphogenesis. **a**, Ventral, dorsal and lateral views of five specimens that are subjectively ordered according to developmental time. Splanchnic mesoderm is represented in blue, and differentiated myocardium is represented in red. The dorsal view shows the progressive dorsal closure of the heart tube and the medial expansion of the SHF. **b**, Limits of the myocardial differentiated domain depicted on the edges of the myocardium midsurface. The green line represents the boundary between the juxta-cardiac field and the myocardium. The blue line represents the boundary between the myocardium and the SHF. **c**, Morphology of the myocardium together with the foregut endoderm and the circulatory system lumen. The different lengths of the aortae seen in the three older specimens are due to variable imaging depth. **d**, Gradual development of the circulatory system in parallel to the foregut endoderm invagination. The ventral view shows the progressive formation of the endocardial lumen, and the lateral view shows the progressive formation of the aortic lumen along the dorsal surface of the foregut endoderm. Abbreviations: JCF, juxta-cardiac field; FGE, foregut endoderm; AP, arterial pole; VP, venous pole; EL, endocardial lumen; AOL, aortic lumen; LSH, left sinus horn; RSH, right sinus horn.</figcaption></figure>

poor classification of the older specimens (Extended Data Fig. 3). The *d1/d2* proportion also provided the most homogeneous distribution of specimens among the 10 groups (Fig. 3d).

The geometries of all the tissues segmented for every specimen classified by stage are shown in Supplementary Figs. 1–8. Groups 1 to 4 represent different stages of CC development, while groups 5 to 8 can be assigned to linear heart tube stages and groups 9 and 10 to heart looping stages. Interestingly, the absolute width of the forming cardiac tube, which approximately coincides with the width of the pericardial cavity, not only does not increase during development but shows a mild trend to reduction as development progresses (Fig. 3e). This aspect indicates that, despite overall embryo growth, the lateral expansion of the pericardial cavity seems restricted. By contrast, the dorsal-ventral and cranio-caudal dimensions of the pericardial cavity and all associated tissues, including the forming heart tube, strongly increase during this period. This aspect correlates well with the extensive deformation that the curve *d1* undergoes,

which is associated with myocardial and splanchnic mesoderm remodeling, in contrast with *d2*, which mostly remains stable.

We also studied the occurrence of developmental events, such as dorsal closure and looping, which are known to progress continuously in live-imaged single specimens¹⁰,¹¹. Heart tube rotation and looping, as described in ref. ¹¹, is appreciable in all specimens of group 10 and most specimens of group 9 (except E42 and E46), without signs of looping in specimens of other groups (Supplementary Figs. 3–6). We also measured the dorsal mesocardial gap as a read out of dorsal closure progression (Extended Data Fig. 4a). This gap is initially wide open in the CC and progressively narrows until disappearing concomitantly with dorsal mesocardium establishment and linear heart tube formation. While dorsal closure appeared predominantly associated with looping and, therefore, was seen in all specimens in group 10+ and most specimens in group 9, we also found some hearts already closed in group 8 (E41) and group 6 (E29), while some specimens at stage 9 were clearly advanced in looping
<figure><img src="image_9.png" /><figcaption>**Fig. 4** | **MSG and morphological variability.** **a,** Illustration of the concept of mean shape for staging group 8 (i) and workflow for obtaining a mean shape and local variability estimation for each staging group (ii). In ii, the left shows the process of surface map computation. The blue shape represents the mesh of the reference specimen of the group. The arrows represent the surface mapping from the reference specimen mesh to all others in the staging group. The right shows the working principle of the surface map computation. The landmarks (lm) identify corresponding locations in the reference and target mesh. Shapes of the group after rigid registration to the reference mesh are also shown (iii). Every vertex in the reference mesh (blue dot) has a corresponding point on the other shapes (red dots). **b,** Calculated myocardium mean shapes for every staging group. The two top rows show the mean shapes from ventral (left) and dorsal (right) views. The two bottom rows show the same mean shapes color coded according to the morphological variability calculated for each vertex. The shape in gray indicates that the low number of specimens at this stage precludes variability calculation. The yellow dotted line and arrows represent the cut trajectory to achieve genus 0 topology in dorsally closed hearts (Methods). Arrowheads indicate the regions of higher variability in medial ventricular regions (solid) and lateral ventricular regions (open). **c,** Representative specimens of each group. Top, ventral and dorsal view of the myocardium shape. Bottom, ventral and dorsal views of the pericardial tissues of the representative embryos.</figcaption></figure>

but quite delayed in dorsal closure (E44 and E48; Supplementary Figs. 3–6). In agreement with these observations, dorsal gap correlation with *d1/d2* was quite noisy (Fig. 3e), indicating that the degree of dorsal closure is not a reliable parameter for staging.

Finally, we wanted to obtain a real-time calibration of the staging system generated. For this, we analyzed 3D+t datasets previously generated in our lab that capture heart tube formation in vivo¹⁰. The limitations of the live-imaging setup did not allow us to measure *d1/d2* or *h/w*; however, we established that the proportion of the ventral protrusion of the primitive chamber to its lateral extension was measurable in live imaging and correlated reasonably well with *d1/d2* (Fig. 3g). Using this approach, we correlated the live-imaging

timing with *d1/d2* (Fig. 3h). This study showed that the available live-imaging data corresponded to the progression from group 1 to group 9, that progression between these stage groups spanned 12 h 42 m and that *d1/d2* and real developmental time showed a linear relationship during this period. These results show that during the live-imaging observation period, and with the caveat of the ex vivo approach, *d1/d2* represents a valid surrogate of developmental time.

**Three-dimensional atlas of heart tube formation with local shape variability.** Given the variability found between specimens in the same staging groups (Supplementary Figs. 1–8), we next aimed to quantify and spatially map the shape variability of the differentiated
myocardium within each of the 10 groups established and to define a consensus geometry for each group (Fig. 4a(i)). For each group, we chose a seed specimen whose *d1/d2* was closest to the average *d1/d2* of the group (Fig. 4a(ii), *y*₁ specimen). We then used a described methodology¹⁸ that uses a set of landmark points to establish a vertex-to-vertex correspondence between all shapes of a group and the seed specimen (Fig. 4a(ii), Extended Data Fig. 4 and Methods). This generated a cluster of corresponding vertices for each position in the mesh. Next, we found the average position of each cluster of corresponding vertices and constructed a mean shape for each group (MSG) composed of the average positions (Fig. 4a(iii)). This generated a collection of 10 MSGs that represents the average evolution of myocardial shape during the formation of the primitive heart tube and initiation of looping (Fig. 4b).

We then computed the per-vertex variability among the specimens of each staging group and plotted it on the MSG collection (Fig. 4a(iii)). In this way, we were able to identify hot spots of variability within each of the 10 different groups (Fig. 4b). In most groups, we found that high variability appeared in the OFTs, inflow tracts (IFTs) and dorsal lips of the myocardium. During primitive ventricular chamber formation, we found that variability concentrated first at the bilateral bulges that start to form the chamber and later at the midline, possibly in relation to the variability in the degree of merging of the two initially bilateral bulges into one. Finally, in group 9, high variability was observed in all structures involved in the looping process (Fig. 4b). Variability was not measured in group 10, given that this group is composed of only two specimens. This collection of MSGs provides information on the standard shapes of different stages of primitive heart tube formation and the local variability of the process.

We also wanted to identify the specimen whose myocardium best represents its stage group, as this would give us access to representative geometries of the other tissues as well. For this, we identified the medoid shape for each group, that is, the geometry with the smallest accumulated difference with respect to all others in the group (Methods). In this way, we identified a concrete specimen as representative shape for each group (RSG; Fig. 4c).

**Early left-right asymmetry of the cardiac IFTs.** Heart tube looping is considered the earliest morphogenetic expression of left-right asymmetry in a mammalian developing organ. Here, we identified obvious signs of heart looping in stages 9 and 10; however, we wanted to explore systematically whether less obvious asymmetries might be present before the looping stage. Studying both the whole specimen collection and the collection of MSGs, we detected a possible recurrent difference in the angle of insertion of the IFTs. To study this aspect in more detail, we defined the direction of insertion of the IFTs via polar coordinates in a local system of orthogonal axes in which the *xy* plane is transversal to the embryo, and the *z* axis is aligned with the embryo AP axis (Fig. 5a, Extended Data Fig. 5 and Methods). The direction of insertion of the IFTs is then defined by two angles: the *θ* angle and the *φ* angle (Fig. 5a). When studying the evolution of these angles

in correlation to the staging parameter *d1/d2*, we found that the *θ* angle progressively diverged between the left and right sides (Fig. 5b,c, ‘Atlas’), while the *φ* angle, on average, remained symmetric (Extended Data Fig. 6a,b). The *θ* angle significantly deviated from symmetry starting at *d1/d2* equal to 1.55, which corresponds to group 6 (Figs. 3f and 5c). The asymmetry of the *θ* angle further increased during the looping stages ($d1/d2 \ge 2$) until reaching a 20° difference (Fig. 5b,c, ‘Atlas’).

Related to the different angle of insertion of the IFTs, we frequently observed that the insertion of the right IFT into the ventricular region forms an acute angle, whereas the similar position on the left side is characterized by a much smoother transition, often even accompanied by a small bulge connecting the ventricle to the IFT (Fig. 5g). In hearts at the looping stage, these asymmetric connections seem to favor the typical rightward ventricular displacement during looping (see the group 9 representative in Fig. 5g), as the right-side acute angle may act as a hinge, favoring rotation toward the right side, while the bulge on the left side may oppose a similar rotation toward the left side.

We concluded that before heart looping is obvious, and starting at linear heart tube stages, left-right morphological asymmetry builds up in the IFT region, manifested in the different angles of insertion of the IFTs.

**Nodal controls early cardiac left-right asymmetry.** To study whether early IFT asymmetry is under the control of the genetic pathway that controls left-right specification in the mammalian embryo, we studied this feature in *Nodal* mutants, which show randomized looping direction¹⁹,²⁰. For this, we studied a collection of 3D images from *Nodal* mutants and their control siblings acquired in the laboratory of S.M.M. First, we determined that the staging systems described here were transferrable between labs, despite the different imaging methods and genetic backgrounds (Methods). We next calculated the *θ* and *φ* angles of the control embryos of this second collection and found progressive asymmetry similar to that observed in the original collection (Fig. 5b). Examining the independent evolution of the left and right *θ* angles in both control embryo collections showed that both tend to diminish with time, but the right one does it at a higher pace than the left one, which produces the progressive asymmetry (Fig. 5d). We then studied conditional *Nodal* mutants (*Hoxb1<sup>Cre/+</sup>; *Nodal<sup>lox/Nud</sup>*) using similar procedures. *Nodal* mutants also showed IFT *θ* angle-specific progressive asymmetry; however, the direction of the asymmetry was reversed with respect to that in control groups (Fig. 5b,c). While in controls, the left-right *θ* angle evolves from zero toward positive values, in *Nodal* mutants, the values evolve from zero toward negative values. This is due to a tendency of *θ* angle to reduce at a higher pace in the left than in the right side in the mutants (Fig. 5e). In all cases, the differences between *Nodal* mutants and the controls start to build up from *d1/d2* values >1.5–1.55. Our results show that early IFT asymmetry is related to heart looping, both in normal development and in mutants, affecting left-right specification (Fig. 5f).

**Fig. 5 | Onset of left-right morphological asymmetry before heart looping.** a, Schemes show the calculated principal directions of the IFT (orange lines) and the polar coordinates system with *θ* and *φ* angles. rsh, right sinus horn; lsh, left sinus horn; vᵢ, vectorial direction of rsh; vᵣ, vectorial direction of lsh. b, Plot showing the evolution of the *θ* angles depending on the *d1/d2* value. Curves represent Padé approximations of every group of specimens. c, Evolution of *P* values as a function of *d1/d2* for the deviation of the *θ* angle from left-right symmetry and for the comparisons between experimental groups, as indicated. A mixed linear model was applied to estimate the asymmetry between both sides at different stages (fixed effects) adjusting by the angle with two-sided *P* value calculation and applying Sidak's correction for multiple tests. d, Separated values of *θ*ᵢ and *θ*ᵣ for the two groups of control specimens. Analysis of covariance has been applied to test whether slopes are significantly different calculating two-tailed *P* values; NS, not significant. e, Separated values of *θ*ᵢ and *θ*ᵣ for the *Nodal*-mutant specimens. The statistical analysis was the same as in d. f, Schematic representation of the different evolution of *θ*ᵢ and *θ*ᵣ angles in wild-type and *Nodal*-mutant embryos. g, Specific wild-type hearts from different groups, highlighting the acute angle between the right IFT insertion and the ventricle (red arrow) and the smooth transition on the opposite side (gray arrow), often accompanied by a small bulge; n, number of specimens showing this feature over the total number of specimens in each staging group. The black arrow illustrates the direction of looping.
**Early left-right asymmetry affects pericardial tissues.** To investigate whether and how the left-right asymmetries detected in the heart tube affect the pericardial tissues, we extended our study to the pericardial cavity mesoderm and the endoderm. To perform this

study, we first determined the medial-sagittal plane of the embryo (Extended Data Fig. 4b) and determined the fitting between the left and right sides of the pericardial mesoderm (Extended Data Fig. 7a) by flipping one onto the other and quantifying the local deformation

a

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /></figure>

b

<figure><img src="image_13.png" /></figure>

c

<figure><img src="image_14.png" /></figure>

d

<figure><img src="image_15.png" /></figure>

f

<figure><img src="image_16.png" /></figure>

e

<figure><img src="image_17.png" /></figure>

g

<figure><img src="image_18.png" /><figcaption>Group 4
n = 2/5</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>Group 5
n = 2/5</figcaption></figure>

<figure><img src="image_20.png" /><figcaption>Group 6
n = 5/7</figcaption></figure>

<figure><img src="image_21.png" /><figcaption>Group 7
n = 5/5</figcaption></figure>

<figure><img src="image_22.png" /><figcaption>Group 8
n = 3/5</figcaption></figure>

<figure><img src="image_23.png" /><figcaption>Group 9
n = 3/5</figcaption></figure>
a

<figure><img src="image_24.png" /></figure>

b

<figure><img src="image_25.png" /></figure>

c

<figure><img src="image_26.png" /></figure>

d

<figure><img src="image_27.png" /></figure>

**Fig. 6** | Dynamic meshes defining the full spatiotemporal domain of heart tube morphogenesis. **a**, Elaboration of a continuous sequence using the MSG shapes. The process starts by mapping the mesh connectivity from the first time point to the next shape. This process is repeated until the connectivity of the first mesh is propagated until the final mesh at $t_{10}$. This process can also be applied using the RSG shapes. **b**, Concept of shape interpolation applied to the three first time points of the RSG shapes. Each row shows the result for the different tissues forming the pericardial cavity. The dots depicted on top of the myocardium shape illustrate the position of a single vertex traveling through the shapes as morphogenesis progresses. **c**, The plot illustrates the estimated trajectory of the vertex highlighted in **b**; myo, myocardium; spl, splanchnic mesoderm; som, somatic mesoderm; par, paraxial mesoderm. **d**, Temporal reference for morphological and physiological events during heart tube formation. Data on the onset of contraction/beating and real-time scale were obtained from in vitro live-imaging videos⁵⁰; HT, heart tube; RV, right ventricle.
required for a perfect match of the geometries (Extended Data Fig. 7b). This study identified asymmetry of the lateral–posterior extensions of the pericardial cavity in which the IFTs are encased. A visualization of the fitting between the right side of the pericardial mesoderm and the flipped left side in the RSG collection suggested a systematic positioning of the left side closer to the midline plane in groups 6+ (Extended Data Fig. 7c). This was confirmed in a morphometric analysis of the whole specimen collection (Extended Data Fig. 7d). This left–right deviation coincides with the direction of the progressive asymmetry of the θ angle detected in the IFTs.

In the endoderm, we observed a deformation of the ventral fore-
gut in correlation with the displacement of the splanchnic meso-
derm layer during dorsal closure (Extended Data Fig. 7e). In early
specimens, this deformation has the form of two bilateral ridges that
coalesce into one, as the two splanchnic mesoderm sides reach the
midline (Extended Data Fig. 7f). From RSG7+, this ridge shows
a deviation from the midline such that its anterior side displaces
rightward, and its posterior side displaces leftward. This displace-
ment is concomitant with a similar tilting of the closing borders of
the splanchnic mesoderm, which produces an asymmetric meso-
cardium, as previously described¹¹. The observations are compatible
with the deformation of the endoderm due to the movement of the
splanchnic mesoderm over its ventral surface.

**A dynamic atlas of the embryonic heart.** To simulate temporal evolution through the stages characterized, we used the 10 MSG shapes generated plus specimen 51 (Supplementary Fig. 5) and the method described in ref. ¹⁸ to generate a smooth temporal transformation between stages (Fig. 6a). To achieve this, we propagated a common mesh connectivity, identifying equivalent vertices from the first to the last shape (Fig. 6a). Once the same mesh connectivity is set on all the MSGs, they were rigidly aligned to the first mesh (t₁). In addition, a second alignment step was taken using the landmarks to register equivalent landmarks across stages. We then applied an interpolation technique (Methods) to estimate the positions of each vertex between consecutive MSGs so that transitions between positions take place smoothly (Fig. 6c) and inserted 30 interpolated shapes between each MSG. This approach generated a dense temporal sequence of 3D meshes of heart tube morphogenesis (Supplementary Video 1). The model generated does not represent the evolution of a concrete specimen but a most probable trajectory of the morphology of the forming heart tube based on averaging the shapes of several hearts at each developmental time.

This approach did not allow us to include the models of the rest of
tissues associated with the forming heart tube because mismatches
are generated at the limits between the MSGs of the different tissues.
To generate a dynamic model that includes all the tissues, we then
used the collection of 10 RSGs (Fig. 6b and Supplementary Videos
2–7). Given that these are specific specimens, their matching to the
surrounding tissues is improved. This model represents a possible
evolution of the morphology of the heart tube and its associated tis-
sues based on different specimens that represent consecutive stages
of heart development. Examples of the interpolated tissues are pro-
vided in Fig. 6b.

Finally, we correlated previously described stages of heart tube
development and several morphological and functional landmarks
of cardiac development to the evolution of *d1/d2* (Fig. 6d) and
incorporated this information to the videos of the four-dimensional
models (Supplementary Videos 2–7).

Discussion

Here, we generated a comprehensive dynamic and realistic repre-
sentation of cardiac tissues during primitive heart tube formation
in a mammalian embryo. In the process of generating a unified
dynamic model using information from a collection of specimens,
we faced the problem of morphological variability during early

cardiogenesis, which represented limitations to describing ‘normal’
development. We thus developed strategies to describe ‘most likely’
geometries for each developmental stage. To achieve this, we first
developed a morphometric staging method independent of subjec-
tive appreciation of specimens, thereby improving the current state
of the art in the field<sup>10,11,14</sup>. Unsupervised clustering and assessment
of parameters that measure the quality of classifications identified
a simple relationship (*d1/d2*) as the best staging parameter. While
we used this parameter to group specimens in relation to their
developmental progression, and this was essential to map variabil-
ity and generation of the 3D+t model, the great advantage of using
*d1/d2* is that it allows the assignment of specimens along a contin-
uum of heart tube morphogenesis. We think this staging method
will be of utility in the precise definition of cardiac developmental
timing when comparing 3D data between different laboratories.
Furthermore, we identified simpler measurements, like the height-
to-width ratio, as very highly correlated with *d1/d2* and accessible to
all labs and 2D–3D imaging approaches.

Ten groups were established according to the variation of *d1/d2*,
which allowed for the statistical analysis of morphology within
groups. The use of a recently developed method¹⁸ for surface com-
parison allowed the definition of mean shapes and the applica-
tion of statistical analyses to each vertex, thereby generating maps
of local variability. We propose that this approach is of high value
when describing mammalian organogenesis because it simultane-
ously provides information about the ‘average’ geometry of a devel-
oping organ and about the regions that are more likely to vary from
specimen to specimen.

Failures in early steps of cardiac development can lead to mis-
alignment of the chambers and great vessels or to different types to
left-right cardiac mispatterning. In particular, we identified strong
variability at the dorsal borders of the myocardium. These dorsal
borders, initially placed laterally, end up fusing at the midline, form-
ing the dorsal mesocardium, which later disappears during looping.
Important mechanical roles have been proposed for the mesocar-
dium during looping¹¹, and the mesocardial myocardium later
localizes to the inner curvature of the looping heart, which plays
essential roles in the alignment of the cardiac chambers and the
great arteries (reviewed in ref. ²¹). The finding of high variability in
this region, therefore, could be related to the onset of malformations
associated with the mesocardial/inner curvature morphogenesis.
Other regions of high variability were the IFTs, which correlate with
the highly dynamic and transient nature of these structures. Finally,
apical regions of the future left ventricle also show high variabil-
ity, possibly associated with asynchrony of the ballooning process²²
between specimens.

Our approach also allowed us to study the left-right symme-
try of the forming heart tube. Heart tube looping is considered as
the earliest morphological evidence of organ left-right patterning
in the embryo. In particular, cardiac left-right patterning defects
are commonly associated with impairing congenital heart disease.
Here, we describe an aspect of left-right asymmetry in the form-
ing heart affecting IFTs and preceding heart looping. The influence
of the IFT region in heart looping has been previously proposed,
mainly through the exertion of asymmetric forces from the left
and right sides of the venous pole, provoked by differences in
cytoskeletal contraction, cell migration and/or proliferation²³⁻²⁵.
We therefore identified the earliest signs of morphological asym-
metry in the mouse embryo localized to the IFTs of the forming
heart. This asymmetry precedes looping and is likely to affect this
process, as any active or reactive force exerted by/on structures
oriented in different angles is predicted to produce asymmetric
forces. In agreement with this view, *Nodal* mutants, which show a
high incidence of heart looping defects, show reversed asymmetry
of the IFT insertion angles. Nonetheless, there is no strict relation-
ship between IFT insertion angles and the looping direction; only
50% of *Nodal* mutants show reversed direction of looping<sup>19,20</sup>, while some of the wild-type specimens show negative left-right θ angles during normal looping. We therefore propose that left-right bias in the orientation of the IFTs depends on the general left-right patterning mechanisms of the embryo and imposes an architectural bias to the forces that determine heart looping direction in interaction with additional determinants. Interestingly, we also found left-right alteration of the pericardial cavity and endoderm, suggesting that the process of heart looping involves a wider asymmetry also affecting the peripheral tissues. Our observations on the relationship between splanchnic mesoderm movements and ventral foregut deformation strongly suggest that the splanchnic mesoderm has an active role in generating left-right asymmetries in the pericardial cavity. In this respect, the described transient expression of *Nodal* in the left splanchnic mesoderm and its influence on cell proliferation and matrix remodeling<sup>20</sup> might represent an important factor in the generation of left-right asymmetry in the pericardial cavity.

The dynamic models generated here represent highly detailed atlases of the different tissues involved in primitive cardiac tube formation. These atlases will be useful for the combination of relevant biological information produced from different experiments or labs into a unique spatial representation. This opens the door for establishing a common space where scientists can share analyses on quantitative aspects representable in the spatiotemporal domain, for example, gene expression patterns. Furthermore, the atlas generated is embedded in a statistical framework that accounts for natural variability, which can be used as a reference for analysis of mutant embryos, either at global or local scales.

Deep understanding of morphogenesis involves identification of the sources and the temporal and spatial distribution of forces that shape developing tissues. Predictive modeling is ultimately required to fully acquire this knowledge<sup>20,26</sup>; however, it critically relies on access to realistic dynamic models of the geometry of the tissues/ organs under study<sup>6</sup>. The dynamic model developed here will be useful for the development of finite-element studies that explain the forces that drive cardiac tube formation, one of the most complex morphogenetic processes in the mammalian embryo, which underlies the high incidence of cardiac malformation in newborns.

## Methods

**Animal model.** The following are mouse alleles used in the study: *Mesp1*<sup>Cre</sup> (ref. 27; MGI:2176467), *Rosa26R*<sup>MrTmG</sup> (ref. 28; MGI:3716464), *ROSA26R*<sup>Tomato</sup> (ref. 29; MGI:6260212), C57BL/6 (Charles River), *Nodal*<sup>box</sup> (ref. 30; MGI:3056345), *Nodal*<sup>mol</sup> (ref. 20) and *Hoxb1*<sup>Cre</sup> (ref. 31; MGI:2668513). Embryos imaged to generate the M.T. collection were on a mixed 129/SvJ × CD1 genetic background. Embryos imaged to generate the S.M.M. collection were on a mixed C57BL/6 × 129/SvJ genetic background. Adult mice were housed in an air-conditioned room with a 12-h light/12-h dark cycle and free access to water and chow diet. A total of 4 males and 10 females were used as progenitors for embryo collection. Animal procedures in the M.T. laboratory were approved by the CNIC Animal Experimentation Ethics Committee, by the Community of Madrid (reference PROEX 144.1/21), and conformed to EU Directive 2010/63EU and Recommendation 2007/526/EC regarding the protection of animals used for experimental and other scientific purposes, enforced in Spanish law under Real Decreto 1201/2005. *Nodal* mutants were housed in the Laboratory of Animal Experimentation and Transgenesis of the SFR Necker, Imagine Campus, Paris. Animal procedures were approved by the ethical committees of the Université de Paris and by the French Ministry of Research.

**Embryo preparation.** Embryos were collected from 00:00 to 08:00 (E8–E8.33). Pregnant mice were killed by CO₂. Embryos were dissected, and the amnion, atlantois and yolk sac were kept, fixed overnight in 2% paraformaldehyde at 4°C, washed three times in PBS and permeabilized in PBS 0.5% Triton X-100 for 30 min. Embryos were incubated overnight at 4°C in PBS with DAPI (1:500) and washed in PBS before embedding in low melting point agarose. Once the agarose polymerized, using a scalpel, a small agarose cube was created (0.5 cm × 0.5 cm × 0.5 cm) containing the embryo. The cube was then immersed in CUBIC I (ref. 32) at 4°C over 3 days, with CUBIC I renovation after the first day. Whole-mount *Nodal* mutants and their controls were fixed at E8.5 in 4% paraformaldehyde and dehydrated in 100% methanol. Embryos were then stained with Hoechst as a nuclear counterstain and transferred to R2 CUBIC clearing reagents.

**Embryo imaging.** Embedded embryos in the laboratory of M.T. were imaged on an inverted microscope through an 80-µm-thick cover glass placed on a perforated culture plate to minimize optical pathway length (Fig. 1a). Given the high sensitivity of CUBIC I to water concentration, a water-saturated atmosphere was created to prevent evaporation. Confocal images were obtained using an SP8 Leica confocal microscope and an HC PL Apo CS2 20 0.75-NA multi-immersion lens corrected to glycerol immersion medium. The following scanning settings were used: xyz mode, bidirectional, 400 Hz; pinhole, 56.7 µm; zoom, 1.15; laser lines used, diode 405 nm (5% power, DAPI signal), 488 nm (power 15%, mGFP signal) and 561 nm (20%, Tomato signal); detectors used, PMT (424 nm–474 nm), HyD (496 nm–529 nm, gain 100% and 583 nm–667 nm, gain 170%). For large specimens, a tile scan was performed, and tiles were stitched using Leica LAS X 3.5.2.18963 software. For the *Nodal*-mutant and control series, 16-bit images were acquired with a Zeiss Z1 light-sheet microscope and a ×20/1.0 objective.

**Image visualization.** Visualization of image stacks was performed in ImageJ 1.53f (ref. 33) and ITK-SNAP 3.8.0.

**Image processing.** Raw images were deconvolved with *Huygens Professional* version 19.10° to compensate for the blurring and noise inherent to the optical system. Deconvolved images were exported in tagged image file format, rectified for axial orientation and saved in NIFTI format. All the image processing steps described hereafter were implemented in Python 3.8 (ref. 34) making use of different libraries for each specific process. The library NiBabel, version 3.1.1 (ref. 35), was used for reading and writing NIFTI images. Images read in this format are represented as array-like structures, defined as an instance of the class ndarray, which is part of the scientific computing package NumPy¹⁶. The library used for processing the segmentation images (Fig. 1f) was scikit-image, version 0.16.2 (ref. 37).

**Image segmentation.** Segmentations were performed using ImageJ in a scaled-down version of the images by applying a diminishing scale factor equal to 2 in each of the three axes. This resulted in lighter files, with xy pixel sizes of 0.76 µm–0.98 µm and a z step of 0.98 µm–4 µm, depending on the specimen. The reduced resolution did not introduce loss of accuracy in the segmentation at the tissue level, because the segmented tissue layers are approximately 20 µm thick. For segmentation of the tissues forming the pericardial cavity and the endocardial lumen, we used the open-source multiplatform ITK-SNAP³⁸. Automatic segmentation based on intensity levels was applied to the Tomato signal, which was present ubiquitously within mesodermal cells. First, we activated the active contour tool and defined the region of interest as the cardiogenic region. Then, we used thresholding as presegmentation mode and adapted the threshold mode and lower and upper values. This creates a new synthetic image called the speed image, which has values close to 1 inside the structure of interest. In the next step, we placed several seeds within the domains of the structure of interest, distributed across the cardiogenic mesoderm. Then, the segmentation was initialized and the seeds extended across the structure of interest. We then performed expert manual delimitation of each of the tissues composing the pericardial cavity and the myocardium and correction of failures detected in the automatic segmentation. For this task, both the mGFP signal, which labels all mesodermal cell membranes, and the Tomato signal, which labels all cell bodies, were used.

The segmentation of the forming aortae and the foregut pocket lumen was done using Mask R-CNN³⁹. The network was trained for the 2D detection of both, the forming aortic lumen and the closed contour defined by the foregut endoderm (as is shown in Fig. 1d). The elaborated dataset was composed of 351 slices taken from random specimens and random planes. In each slice, the contour of every different instance was annotated using ImageJ's ROI manager. The dataset was divided 80/20 for training and validation, respectively. The training and detection were conceived for detection in slices from the frontal view along the ventral-dorsal axis. The performance achieved for the detection of the aortae and foregut pocket was 88% and 91%, respectively. The foregut endoderm was derived from the segmentation of the foregut pocket lumen, performed by selecting negative levels for every channel (Fig. 1d). First, the segmentation image was dilated using the library Scipy, version 1.5.3, using a multidimensional binary dilation. This resulted in an expanded segmented image that covers the surrounding foregut endoderm. Finally, the original segmentation of the foregut pocket was subtracted from the dilated segmentation, resulting in a new segmentation image only representing the foregut endoderm. To determine the limit of the foregut endoderm extension in a reproducible manner across all the specimens, we defined a sphere, centered at the anterior tip of the foregut pocket, with a radius equal to the distance to landmark lm12 (Extended Data Fig. 4c) and calculated the intersection with the foregut endoderm mesh.

**Three-dimensional reconstruction.** After the 2D segmentations were complete, we used the marching cubes algorithm with the library PyMCubes, version 0.1.2, to obtain 3D reconstructions. The resulting meshes were then saved as polygon file format. The meshes were handled in Python by making use of the library Trimesh, version 3.9.19, which was used to filter the meshes, thus creating smooth surfaces. As part of the postprocessing steps, meshes for each of the tissues segmented were first simplified using a simplify quadratic decimation algorithm provided by the
library Open3D, version 0.7.0 (ref. ⁴⁰), to a number of faces equal to 10,000. Then, meshes were checked to make sure that they are closed surface of genus 0 by using Trimesh 3.9.19 and visual inspection using Meshlab 2020.12. The reconstructed surfaces were flipped with respect to the *y* axis of the image (Fig. 1b), the reason being that acquisition in the confocal microscope generates left-right inverted images. By doing this, all surface meshes are oriented as indicated in the system of reference shown in Fig. 1a.

**Embryo collection.** We imaged a total of 52 specimens. All specimens were temporally ordered according to *d1/d2* except two specimens (E51 and E52) that were manually assigned because they represented outliers in the classification (Supplementary Figs. 1–8) and one specimen (E3) in which the partial differentiation of the CC did not allow for the calculation of *d1/d2* and was manually assigned to group 1 by similarity.

Within the collection, there are some exceptions in which we were not able to obtain all the tissues and shapes of interest (Fig. 1a). We list here the name of the embryo, its corresponding staging group and what makes the dataset exceptional.

* E21, group 5: unable to reconstruct the somatic or paraxial mesoderm of this embryo.

* E22, group 5: unable to image the most dorsal part of the foregut endoderm and the aortic lumen in that region.

* E26, group 6: unable to fully reconstruct the foregut endoderm.

* E37, group 8: only myocardium and the endocardial lumen could be reconstructed.

* E41, group 8: unable to segment the splanchnic and somatic mesoderm.

* E51, outlier: unable to segment the foregut endoderm and the aorta.

**Midsurface extraction.** The meshes were first isotropically remeshed¹¹ to a 2-µm edge length and smoothed⁴². Then, the skeleton of the mesh was derived using ref. ¹⁵ and smoothed to eliminate spikes and/or high-frequency details. The former step was repeated until the algorithm created a monolayer mesh. The resulting mesh was then resampled using a Poisson disc sampling and reconstructed using the Ball pivoting algorithm¹¹. After this, the mesh was smoothed and postprocessed to eliminate imperfections (holes, non-manifold edges, non-manifold vertices and so on). The final state of the midsurface comes given by a last isotropic explicit remeshing (2 µm).

**Calculation of distances on the surface of the tissues.** In general, distances between a pair of landmarks were calculated as the geodesic distance between the vertices representing each landmark. The library used for these computations was Networkx, version 2.4 (https://networkx.org). In the cases where the distance between landmarks was defined along the edge of the shape (for example, *d1* and *d2*; Fig. 3a), the distance was calculated by first projecting the landmarks on the midsurface. Then, for each of these points, its closest vertex to the outline (Fig. 3a) was determined. Finally, the geodesic distance was calculated as the accumulated sum of the length of every edge defining the outline between both points.

GraphPad Prism 8.0.1 was used for data representation and statistical analysis.

**Calculation of the dorsal closure gap.** Landmarks 9 and 7 and 8 and 10 were used to define the edges of the dorsal lips of the myocardium (labeled as yellow in Extended Data Fig. 4a). The dorsal gap, *l*, measures the length of the shortest line that joins both lips. The calculation of *l* was achieved by the subdivision of the edge of the lips, generating *N* points in the left side and *K* points in the right side. Then, for every point *j* (left lip) to any point on the right lip, the minimum distance was calculated, *d(p<sub>r</sub>, q)*. Finally, the minimum of the calculated distances was taken as the dorsal gap.

**Surface map computation.** We computed surface maps using the method described in ref. ¹⁸. This method takes two topologically equivalent surface meshes and a sparse set of corresponding points (landmarks) as input. As output, it produces a continuous and bijective (1:1) map of correspondence between the vertices of the two surfaces. The algorithm minimizes intrinsic mapping distortion, which promotes the alignment of geometrically similar surface regions⁴³. Maps are initialized through the layout embedding method presented in ref. ⁴⁴. To establish correspondence between different tissue meshes, we established a sparse set of corresponding points (landmarks) as inputs (Fig. 4a(i) and Extended Data Fig. 4c-f). We provide a software package that implements these procedures, including the vertex-to-vertex mapping tool (SurfaceMapComputation) and a visualization tool (ViewMap) that displays a rigid registration of the shapes. This allows for the adjustment of the landmarks used, the addition of new landmarks and/or trying different initialization algorithms for finding the best mapping between two shapes. The software takes as input parameters the two shapes to be compared, a minimum of four landmarks, number of iterations and a choice of algorithms for the initialization of the map. The initialization algorithm can be automatically selected by the software or manually selected from those described in refs. ⁴⁴,⁴⁵.

**Shape averaging.** We selected one reference shape per group (Fig. 4a(i)). We then computed surface maps¹⁸ from the reference to all other shapes of the group. We used each map to transfer the vertices of the reference mesh to the respective target

surface, establishing a vertex-to-point correspondence. Using this correspondence, we rigidly aligned all shapes of the group to the reference shape to compensate rotational and translational misalignments. We then computed the average 3D position per vertex of the reference mesh, that is, the mean of the reference vertex itself and its corresponding points on the other rigidly aligned shapes (Fig. 4a(ii)). As a measure of local shape variability within the group, we report the standard deviation per set of corresponding vertices given by the square root of the largest eigen value of the covariance matrix (Fig. 4a(ii)).

**Identification of medoid specimens.** The medoid shape of a group is represented by the geometry whose accumulated differences with respect to all others in the group is the smallest (Extended Data Fig. 8, equation 1). First, all the shapes in the group were uniformly scaled to minimize the difference in global size so as to prioritize local morphological variation versus size differences. Then, the accumulated difference *q* for every shape was calculated. This was done by adding all the morphological distances, *D*, between each shape and all the others in the group (Extended Data Fig. 8, equation 2). The morphological distance, *D*, between a pair of shapes was calculated as the average of all the distances between the corresponding vertices of the shapes, *d* (Extended Data Fig. 8, equation 3).

**Sequential mapping to generate the 3D+t atlas.** We established vertex correspondence across the temporal sequence of shapes by computing a surface map between each consecutive pair of shapes. We started from the mesh of the earliest shape and successively connected the vertices to all meshes of the sequence (Fig. 6a). We then rigidly aligned all shapes of the sequence to avoid rotational and translational misalignments between time points. To achieve a high-quality surface approximation between time points, we used a remeshing technique based on ref. ⁴⁶ to obtain a common triangulation that is suited for all keyframe shapes. We computed a target edge length field on each shape of the sequence¹⁷. We then used sequential surface maps to combine all fields via a vertex-wise minimum. In this way, if a shape in the sequence requires a high mesh resolution in some area, this demand will be translated to the whole collection in the common triangulation. We then iteratively improved the common triangulation by a series of edge splits, collapses, flips and tangential Laplacian smoothing⁴⁶ and simultaneously considered the reference mesh embedded on all shapes of the sequence⁴⁸. Finally, to approximate a continuous morphological evolution, we uniformly interpolated 30 new time points between each pair of consecutive shapes in the sequence. We initially determined the intermediate shapes via linear interpolation and then applied a polynomial smoothing filter to elude sharp transitions⁴⁹. We provide a collection of 11 time points plus 30 interpolated time points between keyframes, which represents a total of 311 temporal shapes.

**Surface map computation between shapes of different topology.** The dorsally closed heart and the splanchnic mesoderm represent cases of torus topology (genus 1), while all others have sphere topology (genus 0). This represents a compatibility problem when computing surface maps across the entire sequence. Therefore, we developed a methodology to turn these shapes into sphere topology surfaces by cutting the closed surfaces in a specific region. In the case of the dorsally closed heart, the border between the two myocardial lips is considered as the cut trajectory (Fig. 4b, yellow dotted line and arrow). In the case of the splanchnic, a cut is made along the anterior-posterior axis in the medial domain anterior to the OFT (Fig. 1g, black dotted line and arrow). For the latter, this procedure was done to all the splanchnic shapes of the collection because they all present genus 1 topology.

We processed each cut by individually triangulating the two emerging holes, yielding a closed surface of genus 0. To fuse any visible gap, we slightly extruded both newly triangulated areas (producing a subtle self-intersection of the geometry) and finally remeshed the affected regions.

**Definition of the IFT insertion angles.** The insertion angles are the result of expressing the direction of the IFTs (*v*) via polar coordinates in the system of reference (Fig. 5a). The *θ* angle (Fig. 5a, middle) is defined within the *xy* plane, which normal direction is approximately equal to that of the dorsal-ventral axis of the embryo. To obtain comparable values of *θ* from both IFTs, we flipped the *x* axis of the reference system when computing angles for the right IFT. The *φ* angle (Fig. 5a, right) is defined as the angle through which *v* swings with respect to the cranio-caudal axis. This angle was calculated as the angle formed between *v* and the *z* axis of the reference system. As this angle is defined, it would oscillate around a value of 90°, so we applied an offset of -90° indicated by $θ_{+90}$.

**Application of the staging method between different laboratories.** The S.M.M. collection of control embryos was independently staged using the *h/w* parameter in the S.M.M. and M.T. labs. While the S.M.M. lab used Imaris 9.5.1 on the raw images, the M.T. lab segmented the images and measured on the segmentations. The correlation obtained between the two measurements was very high, indicating that different observers obtained similar measurements from differently treated images (Extended Data Fig. 9b). Next, the M.T. lab calculated *d1/d2* for each of the specimens (Extended Data Fig. 9c). Independently, the specimen segmentations were reconstructed and compared to each of the MSG
collection via surface map computation (as in Extended Data Fig. 8b), which provided a position for each specimen along the MSG sequence with an estimated *d1/d2* value (Extended Data Fig. 9d). For *d1/d2* estimation, the morphometric distances were calculated between the specimen and the MSG of the group with a closest *d1/d2* value plus the two earlier and later MSGs in the temporal sequence (Extended Data Fig. 9e, part 2). The representation of the inverse of these distances (dependent variable) by the *d1/d2* values of the mean shapes (independent variable) was fitted to a Gaussian curve (Extended Data Fig. 9e, part 2). The *d1/d2* value corresponding to the peak fitted curve was identified as the predicted value for the specimen. In this way, we compared for each specimen of the control S.M.M. collection the directly measured *d1/d2* parameter to the *d1/d2* parameter estimated from the matching of the specimen morphology along the MSG sequence. This comparison resulted in a very high linear correlation between the measured and the predicted *d1/d2* values (Extended Data Fig. 9e). This shows that the *d1/d2* parameter and its close surrogate *h/w* can be used to stage collections of embryos of different genetic backgrounds and can be imaged/ measured by different methods and laboratories.

**Renderings.** All the renderings (figures and movies) have been processed in Blender 2.92.

**Reporting Summary.** Further information on research design is available in the Nature Research Reporting Summary linked to this article.

## Data availability

Datasets with the original segmentations, the processed 3D models and the 3D+t models are available at 'Mendeley Data' at https://data.mendeley.com/datasets/t828xhg66k/1. We recommend opening with the open-source software 'Paraview' (www.paraview.org). Instructions for visualizing the 3D images and videos in Paraview with custom selection of tissues, coloring and viewpoint are provided with these datasets.

Datasets with the source confocal images and original segmentations are available at the IDR repository at https://doi.org/10.17867/10000174.

## Code availability

The software package including SurfaceMapComputation and ViewMap is available at Zenodo at https://doi.org/10.5281/zenodo.6390818.

Received: 15 October 2021; Accepted: 6 April 2022;
Published online: 16 May 2022

## References

1. Kelly, R. G., Buckingham, M. E. & Moorman, A. F. Heart fields and cardiac morphogenesis. *Cold Spring Harb. Perspect. Med.* **4**, a015750 (2014).

2. de Bakker, B. S. et al. An interactive three-dimensional digital atlas and quantitative database of human development. *Science* **354**, aag0053 (2016).

3. de Boer, B. A., van den Berg, G., de Boer, P. A., Moorman, A. F. & Ruijter, J. M. Growth of the developing mouse heart: an interactive qualitative and quantitative 3D atlas. *Dev. Biol.* **368**, 203–213 (2012).

4. Faber, J. W., Hagoort, J., Moorman, A. F. M., Christoffels, V. M. & Jensen, B. Quantified growth of the human embryonic heart. *Biol. Open* **10**, bio057059 (2021).

5. Soufan, A. T. et al. Regionalized sequence of myocardial cell growth and proliferation characterizes early chamber formation. *Circ. Res.* **99**, 545–552 (2006).

6. Kawahira, N., Ohtsuka, D., Kida, N., Hironaka, K. I. & Morishita, Y. Quantitative analysis of 3D tissue deformation reveals key cellular mechanism associated with initial heart looping. *Cell Rep.* **30**, 3889–3903 (2020).

7. Mohun, T. J. & Anderson, R. H. 3D anatomy of the developing heart: understanding ventricular septation. *Cold Spring Harb. Perspect. Biol.* **12**, a037465 (2020).

8. Lopez, A. L. III, Wang, S. & Larina, I. V. Embryonic mouse cardiodynamic OCT imaging. *J. Cardiovasc. Dev. Dis.* **7**, 42 (2020).

9. Yue, Y. et al. Long-term, in toto live imaging of cardiomyocyte behaviour during mouse ventricle chamber formation at single-cell resolution. *Nat. Cell Biol.* **22**, 332–340 (2020).

10. Ivanovitch, K., Temiño, S. & Torres, M. Live imaging of heart tube development in mouse reveals alternating phases of cardiac differentiation and morphogenesis. *eLife* **6**, e30668 (2017).

11. Le Garrec, J. F. et al. A predictive model of asymmetric morphogenesis from 3D reconstructions of mouse heart looping dynamics. *eLife* **6**, e28951 (2017).

12. Meilhac, S. M. & Buckingham, M. E. The deployment of cell lineages that form the mammalian heart. *Nat. Rev. Cardiol.* **15**, 705–724 (2018).

13. Susaki, E. A. et al. Whole-brain imaging with single-cell resolution using chemical cocktails and computational analysis. *Cell* **157**, 726–739 (2014).

14. Tyser, R. C. V. et al. Characterization of a common progenitor pool of the epicardium and myocardium. *Science* **371**, eabb2986 (2021).

15. Yoshizawa, S., Belyaev, A. G. & Seidel, H. -P. Free-form skeleton-driven mesh deformations. In *Proc. 8th ACM Symposium on Solid Modeling and Applications* 247–253 (Association for Computing Machinery, 2003).

16. Pedregosa, F. et al. Scikit-learn: machine learning in Python. *J. Mach. Learn. Res.* **12**, 2825–2830 (2011).

17. Rousseeuw, P. J. Silhouettes: a graphical aid to the interpretation and validation of cluster analysis. *J. Comput. Appl. Math.* **20**, 53–65 (1987).

18. Schmidt, P., Campen, M., Born, J. & Kobbelt, L. Inter-surface maps via constant-curvature metrics. *ACM Trans. Graph.* **39**, 119 (2020).

19. Lowe, L. A., Yamada, S. & Kuehn, M. R. Genetic dissection of nodal function in patterning the mouse embryo. *Development* **128**, 1831–1843 (2001).

20. Desgrange, A., Le Garrec, J. F., Bernheim, S., Bonnelykke, T. H. & Meilhac, S. M. Transient Nodal signaling in left precursors coordinates opposed asymmetries shaping the heart loop. *Dev. Cell* **55**, 413–431 (2020).

21. Gittenberger-de Groot, A. C., Bartelings, M. M., Deruiter, M. C. & Poelmann, R. E. Basics of cardiac development for the understanding of congenital heart malformations. *Pediatr. Res.* **57**, 169–176 (2005).

22. Christoffels, V. M. et al. Chamber formation and morphogenesis in the developing mammalian heart. *Dev. Biol.* **223**, 266–278 (2000).

23. Ocaña, O. H. et al. A right-handed signalling pathway drives heart looping in vertebrates. *Nature* **549**, 86–90 (2017).

24. Kidokoro, H., Okabe, M. & Tamura, K. Time-lapse analysis reveals local asymmetrical changes in C-looping heart tube. *Dev. Dyn.* **237**, 3545–3556 (2008).

25. Voronov, D. A., Alford, P. W., Xu, G. & Taber, L. A. The role of mechanical forces in dextral rotation during cardiac looping in the chick embryo. *Dev. Biol.* **272**, 339–350 (2004).

26. Sharpe, J. Computer modeling in developmental biology: growing today, essential tomorrow. *Development* **144**, 4214–4225 (2017).

27. Saga, Y. et al. MesP1 is expressed in the heart precursor cells and required for the formation of a single heart tube. *Development* **126**, 3437–3447 (1999).

28. Muzumdar, M. D., Tasic, B., Miyamichi, K., Li, L. & Luo, L. A global double-fluorescent Cre reporter mouse. *Genesis* **45**, 593–605 (2007).

29. Madisen, L. et al. Transgenic mice for intersectional targeting of neural sensors and effectors with high specificity and performance. *Neuron* **85**, 942–958 (2015).

30. Lu, C. C. & Robertson, E. J. Multiple roles for Nodal in the epiblast of the mouse embryo in the establishment of anterior–posterior patterning. *Dev. Biol.* **273**, 149–159 (2004).

31. Arenkiel, B. R., Gaufo, G. O. & Capecchi, M. R. Hoxb1 neural crest preferentially form glia of the PNS. *Dev. Dyn.* **227**, 379–386 (2003).

32. Matsumoto, K. et al. Advanced CUBIC tissue clearing for whole-organ cell profiling. *Nat. Protoc.* **14**, 3506–3537 (2019).

33. Schindelin, J. et al. Fiji: an open-source platform for biological-image analysis. *Nat. Methods* **9**, 676–682 (2012).

34. Van Rossum, G. D. & Fred L. Python 3 Reference Manual (CreateSpace, 2009).

35. Brett, M. et al. nipy/nibelabel: 3.1.1. Zenodo https://doi.org/10.5281/ZENODO.3924343 (2020).

36. Harris, C. R. et al. Array programming with NumPy. *Nature* **585**, 357–362 (2020).

37. van der Walt, S. et al. Scikit-image: image processing in Python. *PeerJ* **2**, e453 (2014).

38. Yushkevich, P. A. et al. User-guided 3D active contour segmentation of anatomical structures: significantly improved efficiency and reliability. *NeuroImage* **31**, 1116–1128 (2006).

39. He, K., Gkioxari, G., Dollár, P. and Girshick, R. Mask R-CNN. 2017 IEEE International Conference on Computer Vision (ICCV) 2980–2988 (2017).

40. Zhou, Q.-Y., Park, J. & Koltun, V. Open3D: a modern library for 3D data processing. Preprint at https://doi.org/10.48550/arXiv.1801.09847 (2018).

41. Cignoni, P. et al. MeshLab: an open-source mesh processing tool. In *Italian Chapter Conference 2008* (Ed. Scarano, V., De Chiara, R. & Erra, U.) (The Eurographics Association, 2008).

42. Taubin, G. Curve and surface smoothing without shrinkage. In *Proc. IEEE International Conference on Computer Vision* pp. 852–857 (IEEE, 1995).

43. Schmidt, P., Born, J., Campen, M. & Kobbelt, L. Distortion-minimizing injective maps between surfaces. *ACM Trans. Graph.* **38**, 156 (2019).

44. Born, J., Schmidt, P. & Kobbelt, L. Layout embedding via combinatorial optimization. *Comput. Graphics Forum* **40**, 277–290 (2021).

45. Schreiner, J., Asirvatham, A., Praun, E. & Hoppe, H. Inter-surface mapping. *ACM Trans. Graph.* **23**, 870–877 (2004).

46. Botsch, M. & Kobbelt, L. A remeshing approach to multiresolution modeling. In *Proc. 2004 Eurographics/ACM SIGGRAPH Symposium on Geometry Processing* pp. 185–192 (Association for Computing Machinery, 2004).

47. Dunyach, M., Vanderhaeghe, D., Barthe, L. & Botsch, M. Adaptive remeshing for real-time mesh deformation. In *Eurographics 2013—Short Papers* (Ed. Otaduy, M. -A. & Sorkine, O.) (The Eurographics Association, 2013).

48. Yang, Y., Zhang, W.-X., Liu, Y., Liu, L. & Fu, X.-M. Error-bounded compatible remeshing. *ACM Trans. Graph.* **39**, 113 (2020).
49. Luo, J., Ying, K. & Bai, J. Savitzky-Golay smoothing and differentiation filter for even number data. *Signal Process.* **85**, 1429–1434 (2005).

50. Ivanovitch, K., Esteban, I. & Torres, M. Growth and morphogenesis during early heart development in amniotes. *J Cardiovasc. Dev. Dis.* **4**, 20 (2017).

## Acknowledgements

We thank members of the Torres laboratory for fruitful discussions and advice on this work. Microscopy was conducted at the Microscopy & Dynamic Imaging, CNIC, ICTS-ReDib, co-funded by MCIN/AEI /10.13039/501100011033 and FEDER Una manera de hacer Europa" (#ICTS-2018-04-CNIC-16). We are grateful to F. Sánchez-Cabo for helpful advice on the statistics. I.E. acknowledges attendance to the Computational Image Analysis (CIAN) course at the MBL, Woods Hole, MA (funded by NIH grant R25 GM103792). We thank S. Bernheim, T. Holm Bønnelykke, L. Guillemot and M. Franco for technical assistance and the Imaging platform and the LEAT animal facility of the SFR Necker. The following provided funding for this study: grant PGC2018-096486-B-100 from the Spanish Ministerio de Ciencia e Innovación to M.T., grant H2020-MSCA-ITN-2016-722427 from the EU Horizon 2020 program to M.T., S2017/BMD-3875 from the Comunidad de Madrid to M.T., predoctoral fellowship BFU2015-71519-P from the Spanish Ministry of Economy, Industry and Competitiveness (MEIC) to I.E., partial funding from the German Research Foundation within the Gottfried Wilhelm Leibniz program the Excellence Strategy of the Federal Government and the Länder to L.K., support from the core funding from the Institut Pasteur to S.M.M. and state fundings from the Agence Nationale de la Recherche under the 'Investissements d'avenir' program (ANR-10-IAHU-01) to S.M.M. S.M.M. is an INSERM research director. The CNIC is supported by the Ministerio de Ciencia e Innovación and the Pro CNIC Foundation and is a Severo Ochoa Center of Excellence (CEX2020-001041-S).

## Author contributions

I.E. and M.T. conceptualized the project. S.T., I.E., P.S. and A.D. developed the methodology. I.E. and M.R. performed formal analyses. I.E. and M.T. wrote the original draft. I.E., P.S., S.M.M., L.K. and M.T. reviewed and edited the final paper. S.M.M., L.K. and M.T. supervised the project. S.M.M., L.K. and M.T. acquired funding.

## Competing interests

The authors declare no competing interests.

## Additional information

**Supplementary information** The online version contains supplementary material available at https://doi.org/10.1038/s44161-022-00065-1.

**Correspondence and requests for materials** should be addressed to Miguel Torres.

**Peer review information** *Nature Cardiovascular Research* thanks Badri Roysam, Gonzalo del Monte Nieto and the other, anonymous, reviewers for their contribution to the peer review of this work.

**Reprints and permissions information** is available at www.nature.com/reprints.

**Publisher's note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_28.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long

as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2022
<figure><img src="image_29.png" /></figure>

<figure><img src="image_30.png" /><figcaption>Parameter 7</figcaption></figure>

<figure><img src="image_31.png" /></figure>

<figure><img src="image_32.png" /><figcaption>Parameter 8</figcaption></figure>

<figure><img src="image_33.png" /></figure>

<figure><img src="image_34.png" /><figcaption>Parameter 9</figcaption></figure>

<figure><img src="image_35.png" /></figure>

<figure><img src="image_36.png" /><figcaption>Parameter 10</figcaption></figure>

<figure><img src="image_37.png" /></figure>

<figure><img src="image_38.png" /><figcaption>Profile</figcaption></figure>

<figure><img src="image_39.png" /></figure>

<figure><img src="image_40.png" /><figcaption>Profile</figcaption></figure>

<figure><img src="image_41.png" /></figure>

<figure><img src="image_42.png" /><figcaption>Parameter 13</figcaption></figure>

<figure><img src="image_43.png" /><figcaption>fl</figcaption></figure>

<figure><img src="image_44.png" /><figcaption>Parameter 15</figcaption></figure>

<figure><img src="image_45.png" /></figure>

<figure><img src="image_46.png" /><figcaption>Parameter 16</figcaption></figure>

<figure><img src="image_47.png" /></figure>

<figure><img src="image_48.png" /><figcaption>Parameter 17</figcaption></figure>

<figure><img src="image_49.png" /></figure>

<figure><img src="image_50.png" /><figcaption>Parameter 18</figcaption></figure>

<figure><img src="image_51.png" /></figure>

<figure><img src="image_52.png" /><figcaption>Parameter 19</figcaption></figure>

<figure><img src="image_53.png" /></figure>

<figure><img src="image_54.png" /><figcaption>Parameter 20</figcaption></figure>

**Extended Data Fig. 1** | Additional parameters explored for staging. Schematic representation of additional parameters explored for the staging system.
a

<figure><img src="image_55.png" /></figure>

b

<figure><img src="image_56.png" /></figure>

c

<figure><img src="image_57.png" /></figure>

d

<figure><img src="image_58.png" /></figure>

e

<figure><img src="image_59.png" /></figure>

f

<figure><img src="image_60.png" /></figure>

**Extended Data Fig. 2** | Silhouette coefficient plot for the different staging parameters explored, classified using k-means and number of groups equal to 10. **a–f**, Bar graphs showing the silhouette coefficient value (s) for each embryo in the collection classified into 10 groups, according to the different staging parameters. The graphs also show the average value (⟨s⟩), which is maximal for parameter 1. Graphs are color-coded according to the table shown in Fig. 3D. Schemes above the graphs show a representation of staging parameter used.
<figure><img src="image_61.png" /></figure>

<figure><img src="image_62.png" /></figure>

<figure><img src="image_63.png" /></figure>

<figure><img src="image_64.png" /></figure>

<figure><img src="image_65.png" /></figure>

**Extended Data Fig. 3** | Correlations between staging parameters unveil tight relationships between the evolution of different morphometric aspects during heart development. **a, b, c, d, e**, Plots show the correlation of staging parameters 2-6 versus staging parameter 1, *d*1/*d*2.
a

<figure><img src="image_66.png" /></figure>

DMG calculation
$DMG = \min(d(p_j, q)) \forall j \in J$

<figure><img src="image_67.png" /></figure>

Distance from $p_j$ to the other lip
$d(p_j, q) = \min(d(p_j, q_k)) \forall k \in K$

d

<figure><img src="image_68.png" /></figure>

<figure><img src="image_69.png" /></figure>

<figure><img src="image_70.png" /></figure>

<figure><img src="image_71.png" /></figure>

b

<figure><img src="image_72.png" /></figure>

c

<figure><img src="image_73.png" /></figure>

f

<figure><img src="image_74.png" /></figure>

e

<figure><img src="image_75.png" /></figure>

Extended Data Fig. 4 | See next page for caption.
**Extended Data Fig. 4** | Measurement of dorsal closure and definition of landmarks and curves for surface map computation and shape registration.

a, Schematic representation of the procedure for estimating the dorsal mesocardial gap (DMG). b, This panel shows the process for estimating the mid-sagittal plane of the embryo. The right part shows the tissues of the pericardial cavity cut in half through the plane. A minimum of 3 mid-sagittal points determined in the raw images were used to calculate the mid-sagittal plane. c, The panel shows the basic set of landmarks defined for the myocardium. Im1, Im2, Im3 and Im4 were set at the lateral ends of the IFTs. Im5 and Im6 are located at the point where the IFTs join the ventricle. Im7, Im8, Im9 and Im10 are located at the ends of the closing dorsal lips. Im11 is located at the cranial-most point of the intersection between the mid-sagittal plane and the OFT. Im12 is located at the intersection between the mid-sagittal plane and the line between the myocardium and the juxta-cardiac field. The curves extracted from the edge of the mid-surface are: curve1, from Im3 to Im4; curve2, from Im1 to Im2; curve3, from Im1 to Im3 (left) and Im2 to Im4 (right).

d, The panel shows the set of landmarks defined in the splanchnic mesoderm. Im3, Im4 and Im12 are the nearest vertex to those in the myocardial landmarks. Im13 is defined at the middle position of the intersection between the splanchnic and the somatic mesoderm. Im15 is at the cranial-most point of the intersection between the mid-sagittal plane and the splanchnic mesoderm. Im14 is obtained by calculating the midpoint of the geodesic curve that joins Im12 and Im13. Lm16 and Im17 are located at the bilateral most caudal points of the splanchnic. Im24 and Im25 are at the midpoint of the geodesic curve that joins Im16 with Im4 and Im17 with Im3, respectively. The curves extracted from the edge of the mid-surface are: curve4, at the border between the splanchnic and somatic mesoderm; curve5, a bilateral curve joining Im24, Im15 and Im25; curve6, joins Im24 and Im25 running equidistant between curve2 and curve4; curve7 is similar to curve5, but running equidistant between curve1 and curve5. e, The panel shows the basic set of landmarks defined in the somatic mesoderm. Im13 and Im15 are located at the midpoint of the outer edge of the somatic mesoderm. Im26 is at the midpoint of the geodesic curve that joins Im13 and Im15. Lm18 and Im19 are the most caudal bilateral points of the somatic mesoderm. The curves extracted from the edge of the mid-surface are: curve8, at the border between the somatic and the splanchnic mesoderm; curve9, from Im18 to Im19, on the border between somatic and the paraxial mesoderm; curve10, starts at Im18 and ends at Im19 running equidistant between curve8 and curve9. f, The panel shows the basic set of landmarks defined in the head paraxial mesoderm. Im20 and Im21 are at the bilateral caudal-most points of the paraxial mesoderm. Im22 and Im23 are defined at the bilateral cranial-most points of the paraxial domain. The curves extracted from the edge of the mid-surface are: curve11, running between the paraxial and somatic mesoderm; curve12, running between the paraxial and the splanchnic mesoderm; curve13, running equidistant between curve11 and curve12.
a

1.

<figure><img src="image_76.png" /></figure>

Using the Mid-surface, Im5 and Im6
are projected to the **boundary**

2.

<figure><img src="image_77.png" /></figure>

The former points are then
projected to the **other boundary**

3.

<figure><img src="image_78.png" /></figure>

The former points are used to
define the limits of the inflow tracts

4.

<figure><img src="image_79.png" /></figure>

Landmarks 1, 2, 3 and 4 define
**inflow segments**

5.

<figure><img src="image_80.png" /></figure>

The former points are used to
define the **ventricle boundary**

6.

<figure><img src="image_81.png" /></figure>

Vertices of the boundary are fitted
to a line, which defines the y axis

7.

<figure><img src="image_82.png" /></figure>

Segments are equidistantly divided
and joined by **lines** from top to bottom

8.

<figure><img src="image_83.png" /></figure>

Average of the normals of each
first **line** gives the x axis

9.

<figure><img src="image_84.png" /></figure>

Cross product of x and y
gives z axis

10.

<figure><img src="image_85.png" /></figure>

Joining the center of the lines
represents **1D inflows tracts**

11.

<figure><img src="image_86.png" /></figure>

Linear fitting of the **1D inflows tracts**
gives the **inflows principal directions**

12.

<figure><img src="image_87.png" /></figure>

Inflow tracts principal directions
and inherent system of reference

**Extended Data Fig. 5** | Proposed algorithm for systematic calculation of a myocardium inherent system of reference and calculation of IFT insertion direction. Steps followed for determining the principal direction of the IFTs and a system of reference based on morphology.
<figure><img src="image_88.png" /></figure>

**Extended Data Fig. 6** | The inflow tract shows $\varphi$ angle left-right symmetry. **a**, Plot showing the evolution of $\varphi$ angles L-R symmetry depending on the $d1/d2$ parameter. Curves represent Padé approximations for each experimental group. **b**, Evolution of p-values as a function of $d1/d2$ for the deviation of the $\theta$ angle from L-R symmetry and for the comparisons between experimental groups, as indicated. A mixed linear model was applied to estimate the asymmetry between both sides at different stages (fixed effects) adjusted by the angle with two-sided p-value calculation.
a

<figure><img src="image_89.png" /></figure>

Pericardial
cavity

<figure><img src="image_90.png" /></figure>

Myocardium

<figure><img src="image_91.png" /></figure>

<figure><img src="image_92.png" /></figure>

Pericardial
mesoderm

b

<figure><img src="image_93.png" /></figure>

Right

Left

<figure><img src="image_94.png" /></figure>

<figure><img src="image_95.png" /></figure>

Left

Right

<figure><img src="image_96.png" /></figure>

c

Left

Right

<figure><img src="image_97.png" /></figure>

d

<figure><img src="image_98.png" /></figure>

e

<figure><img src="image_99.png" /></figure>

<figure><img src="image_100.png" /></figure>

f

<figure><img src="image_101.png" /></figure>

<figure><img src="image_102.png" /></figure>

Extended Data Fig. 7 | See next page for caption.
**Extended Data Fig. 7** | Signs of early asymmetry in the pericardial mesoderm. **a**, Reconstruction of the pericardial mesoderm. **b**, Maps of local deformation distribution resulting from mapping left and right sides of the pericardial mesoderm to each other. Specimens shown are the RSGs from ventral and dorsal views. Dotted magenta line indicates the midplane intersection. The color coding corresponds to the relative levels of morphological variation within each specimen. **c**, Ventral view of the pericardial mesoderm of RSGs (orange), including the view of the left side flipped with respect to the midplane onto the right side. **d**, Plot shows the characterization of the misalignment between the two sides of the pericardial mesoderm. A mixed linear model was applied to estimate the asymmetry between both sides at different stages (fixed effects) adjusted by the *l* distance with two-sided *p*-value calculation. The schemes on the right show how the *l* distance was determined. **e**, Ventral view of the foregut pocket endoderm in the RSGs. **f**, A similar representation of the foreguts in *e*, indicating the midplane intersection (pink dotted lines), the ridges detected, first bilateral and later fused medially (blue), and the *d1* line, at the border between the myocardium and the splanchnic mesoderm (red). **n**, normal direction of the midplane.
Group shapes: $y_1 \ y_2 \ y_3 \ y_4 \boxed{y_5} \rightarrow y_{medoid} = \arg \min_{x_i \in \chi} q_{x_i} \quad (1)$

Medoid shape calculation

1. Uniform scaling

<figure><img src="image_103.png" /></figure>

2. Calculation of $q$:

<figure><img src="image_104.png" /></figure>

3. Distance between two shapes, $D$:

<figure><img src="image_105.png" /></figure>

$$D_{x_1,x_5} = \frac{1}{|V|} \sum_{v \in V} d_{x_1,v, x_5,v} \quad (3)$$

<figure><img src="image_106.png" /></figure>

**Extended Data Fig. 8 | Details of the calculation of the medoid shape of a group.** Figure shows the steps to determine the medoid shape of a group, following the terms included in equation 1. Shapes named and enumerated in the figure correspond to group 8. **y**, indicates an individual shape of the group; **Y**, indicates the set of shapes in the group; **x**<sub>i</sub>, indicates a shape of the group, after uniform scaling; **X**, indicates the set of shapes in the group after uniform scaling; **q**, indicates the accumulated sum of morphological distances between a shape and all the other in the group; **n**, number of shapes in the group; **D**, morphological distance between a pair of shapes; **v**, index of a specific vertex; **V**, set of vertices.
<figure><img src="image_107.png" /><figcaption>b</figcaption></figure>

<figure><img src="image_108.png" /><figcaption>c</figcaption></figure>

<figure><img src="image_109.png" /><figcaption>d</figcaption></figure>

<figure><img src="image_110.png" /><figcaption>e</figcaption></figure>

<figure><img src="image_111.png" /></figure>

**Extended Data Fig. 9** | Pipeline for staging the specimens from the study of *Nodal* loss of function. **a**, Steps taken for the image acquisition and segmentation of the *Nodal* mutant and control specimen series. **b**, Correlation between the *h/w* staging parameter independently measured in the SMM and MT labs. **c**, Mapping of the *d*1/*d*2 values measured in SMM's control specimen collection. **d**, Procedure for the morphological mapping of the SMM's control specimen collection to the atlas MSG collection. **e**, Correlation between the experimentally determined *d*1/*d*2 parameter and the predicted *d*1/*d*2 obtained from the morphological mapping of the SMM's control specimen collection to the atlas MSG collection.
# Reporting Summary

Nature Portfolio wishes to improve the reproducibility of the work that we publish. This form provides structure for consistency and transparency in reporting. For further information on Nature Portfolio policies, see our Editorial Policies and the Editorial Policy Checklist.

## Statistics

For all statistical analyses, confirm that the following items are present in the figure legend, table legend, main text, or Methods section.

n/a

☐ The exact sample size (n) for each experimental group/condition, given as a discrete number and unit of measurement

☐ A statement on whether measurements were taken from distinct samples or whether the same sample was measured repeatedly

☐ The statistical test(s) used AND whether they are one- or two-sided
Only common tests should be described solely by name; describe more complex techniques in the Methods section.

☑ A description of all covariates tested

☐ A description of any assumptions or corrections, such as tests of normality and adjustment for multiple comparisons

☐ A full description of the statistical parameters including central tendency (e.g. means) or other basic estimates (e.g. regression coefficient) AND variation (e.g. standard deviation) or associated estimates of uncertainty (e.g. confidence intervals)

☐ For null hypothesis testing, the test statistic (e.g. F, t, r) with confidence intervals, effect sizes, degrees of freedom and *P* value noted
Give *P* values as exact values whenever suitable.

☑ For Bayesian analysis, information on the choice of priors and Markov chain Monte Carlo settings

☑ For hierarchical and complex designs, identification of the appropriate level for tests and full reporting of outcomes

☑ Estimates of effect sizes (e.g. Cohen's *d*, Pearson's *r*), indicating how they were calculated

Our web collection on statistics for biologists contains articles on many of the points above.

## Software and code

Policy information about availability of computer code

**Data collection** Leica LAS X 3.5.2.18963 software was used to stitch tile-scan stacks. Huygens Professional version 19.10 was used to deconvolve raw images.

**Data analysis** Graphpad Prism 8.0.1 was used for data representation and statistical analysis. ITK-SNAP 3.8.0 was used for image segmentation. ImageJ 1.53f was used for image visualization and processing. Paraview 5.4 was used for mesh visualization. NiBabel 3.1.1 was used for reading and writing Nifti image format (segmentation images). All the image processing steps have been implemented in Python 3.8. Scikit-image 0.16.2 was used for denoising and filtering images. Trimesh 3.9.19 and Meshlab 2020.12 were used for mesh analysis and processing. PyMCubes 0.1.2 was used to reconstruct in 3D the segmentation images. The library "Manifold Approximation of 3D Medial Axis" (Shin Yoshizawa, 2011 release) was used to skeletonize the meshes and produce the mid-surfaces. Imaris 9.5.1 has been used by S.M.M. to measure the h/w ratios. A new software package for Intersurface maps computation is provided through the open repository Zenodo.

For manuscripts utilizing custom algorithms or software that are central to the research but not yet described in published literature, software must be made available to editors and reviewers. We strongly encourage code deposition in a community repository (e.g. GitHub). See the Nature Portfolio guidelines for submitting code & software for further information.

## Data

Policy information about availability of data

All manuscripts must include a data availability statement. This statement should provide the following information, where applicable:

- Accession codes, unique identifiers, or web links for publicly available datasets

- A description of any restrictions on data availability

- For clinical datasets or third party data, please ensure that the statement adheres to our policy

DATA AVAILABILITY
Datasets with the original segmentations, the processed 3D models and the 3D+t models are available at Mendeley Data: https://data.mendeley.com/datasets/t828xhg66k/1 (temporary access at: https://data.mendeley.com/datasets/t828xhg66k/draft?aa=33d4dd5d-d49c-4929-b2d5-b41b37d5d418). We recommend opening with the open-source software "Paraview" (www.paraview.org). Instructions for visualizing the 3D images and Videos in ParaView with custom selection of tissues, colouring and viewpoint are provided with these datasets.

Datasets with the source confocal images and original segmentations are available at the IDR repository: https://doi.org/10.17867/10000174

## CODE AVAILABILITY

The software package including SurfaceMapComputation and ViewMap is available at Zenodo repository: https://doi.org/10.5281/zenodo.6390818

# Field-specific reporting

Please select the one below that is the best fit for your research. If you are not sure, read the appropriate sections before making your selection.

<figure><img src="image_112.png" /></figure>

Life sciences

<figure><img src="image_113.png" /></figure>

Behavioural & social sciences

<figure><img src="image_114.png" /></figure>

Ecological, evolutionary & environmental sciences

For a reference copy of the document with all sections, see [nature.com/documents/nr-reporting-summary-flat.pdf](https://nature.com/documents/nr-reporting-summary-flat.pdf)

# Life sciences study design

All studies must disclose on these points even when the disclosure is negative.

Sample size

The main sample collection (52 specimens) was used to generate a sequence of developmental stages. The aim was to define a specimen collection with sufficient temporal density to describe the morphological evolution of tissues. In this sense, the determination of the number of specimens is arbitrary. We chose a temporal density within the range of previous developmental Biology studies aiming to generate similar data for other developing organs. The study that has used a temporally denser collection of specimens for mouse development is that described in Musy et al. (Development (2018) 145 (7): dev154856). In that study, they collected embryos every 8 minutes on average. Their study then concluded that it was not possible to detect developmental differences below 2h of difference between embryos. We therefore set our temporal density to 1 specimen every 20 minutes, which would provide around 6 embryos each 2 hours. This was determined as the best compromise between number of specimens and statistical description of morphological trajectories of the developing heart

Data exclusions

No data were excluded

Replication

A second collection was analyzed containing 39 specimens, including controls and mutants. The sample size was estimated from the data obtained with the main collection. This second collection (replication) validated the staging system defined by the first collection. No other hypotheses were tested and only this replication was performed.

Randomization

We did not compare groups of specimens under different treatments or conditions, so randomization is not applicable to any of our studies.

Blinding

Groups were defined by genotype and data acquisition was not blinded regarding group allocation. Nonetheless, for comparison between groups, data were extracted automatically from the 3D images, so no subjective judgment was involved at any step

# Reporting for specific materials, systems and methods

We require information from authors about some types of materials, experimental systems and methods used in many studies. Here, indicate whether each material, system or method listed is relevant to your study. If you are not sure if a list item applies to your research, read the appropriate section before selecting a response.

## Materials & experimental systems

<table><tr><td>n/a</td><td>Involved in the study</td></tr><tr><td>☑</td><td>Antibodies</td></tr><tr><td>☑</td><td>Eukaryotic cell lines</td></tr><tr><td>☑</td><td>Palaeontology and archaeology</td></tr><tr><td>☐</td><td>Animals and other organisms</td></tr><tr><td>☐</td><td>Human research participants</td></tr><tr><td>☐</td><td>Clinical data</td></tr><tr><td>☐</td><td>Dual use research of concern</td></tr></table>

## Methods

<table><tr><td>n/a</td><td>Involved in the study</td></tr><tr><td>☑</td><td>ChIP-seq</td></tr><tr><td>☑</td><td>Flow cytometry</td></tr><tr><td>☑</td><td>MRI-based neuroimaging</td></tr></table>

## Animals and other organisms

Policy information about studies involving animals; ARRIVE guidelines recommended for reporting animal research

Laboratory animals

A total of 91 laboratory mouse (Mus musculus) embryos between days 7 and 8.5 of gestation were used. This corresponded to 52 embryos with fluorescent reporter transgenes and Cre lines: Mesp1-Cre, ROSA26R-mTmG and ROSA26R-Tomato in the M.T. laboratory and 39 unlabeled embryos in the S.M.M. laboratory. A total of 4 males and 10 females females were used as progenitors
for embryo collection.

Mouse alleles used in the manuscript are listed including bibliographic references and allele identities at the 'Mouse Genome Informatics' database (MGI, http://www.informatics.jax.org/). Mesp1cre ('MGI:2176467'), Rosa26RmTmG (MGI:3716464), ROSA26RTomato (MGI:6260212), C57BL/6 (Charles River), Nodalflox (MGI:3056345), Nodalnul, Hoxb1Cre (MGI:2668513). Adult mice were housed in an air-conditioned room with a 12h light/dark cycle and free access to water and chow diet. Embryos imaged to generate the MT collection were maintained on a mixed 129/SvJ x CD1 genetic background. Embryos imaged to generate the SMM collection were maintained on a mixed C57BL/6 x 129/SvJ genetic background. All animal procedures in the MT laboratory were approved by the CNIC Animal Experimentation Ethics Committee, by the Community of Madrid (Ref. PROEX 144.1/21) and conformed to EU Directive 2010/63EU and Recommendation 2007/526/EC regarding the protection of animals used for experimental and other scientific purposes, enforced in Spanish law under Real Decreto 1201/2005. Nodal mutants were housed in the Laboratory of Animal Experimentation and Transgenesis of the SFR Necker, Imagine Campus, Paris. Animal procedures were approved by the ethical committees of the Université de Paris and by the French Ministry of Research.

Wild animals

The study does not include wild animals

Field-collected samples

The study does not include field-collected samples

Ethics oversight

All animal procedures in the MT laboratory were approved by the CNIC Animal Experimentation Ethics Committee and by the Community of Madrid (Ref. PROEX 144.1/21) and conformed to EU Directive 2010/63EU and Recommendation 2007/526/EC regarding the protection of animals used for experimental and other scientific purposes, enforced in Spanish law under Real Decreto 1201/2005. Nodal mutants were housed in the Laboratory of Animal Experimentation and Transgenesis of the SFR Necker, Imagine Campus, Paris, and animal procedures were approved by the ethical committees of the Université de Paris and by the French Ministry of Research.

Note that full information on the approval of the study protocol must also be provided in the manuscript.