Quantification of gene expression patterns to reveal the origins of abnormal morphogenesis

Abstract

The earliest developmental origins of dysmorphologies are poorly understood in many congenital diseases. They often remain elusive because the first signs of genetic misregulation may initiate as subtle changes in gene expression, which are hard to detect and can be obscured later in development by secondary effects. Here, we develop a method to trace back the origins of phenotypic abnormalities by accurately quantifying the 3D spatial distribution of gene expression domains in developing organs. By applying Geometric Morphometrics to 3D gene expression data obtained by Optical Projection Tomography, we determined that our approach is sensitive enough to find regulatory abnormalities that have never been detected previously. We identified subtle but significant differences in the gene expression of a downstream target of a Fgfr2 mutation associated with Apert syndrome, demonstrating that these mouse models can further our understanding of limb defects in the human condition. Our method can be applied to different organ systems and models to investigate the etiology of malformations.

Keywords: Apert syndrome; Geometric Morphometrics (GM); Optical Projection Tomography (OPT); developmental biology; developmental defects; limb development; mouse; whole‐mount‐in‐situ hybridization (WMISH).

Conflict of interest statement

NM, RM, JS, SM, MY, AR, JS, LR, KK, JR, JS No competing interests declared

© 2018, Martínez-Abadías et al.

Figures

Figure 1.. Quantitative size and shape comparison of forelimb bones in Fgfr2+/P253R newborn mice (P0) and unaffected littermates.

(A) Mouse skeleton at P0. 3D isosurface reconstruction of the skeleton of an unaffected littermate obtained from a high-resolution µCT scan. (B–F) Anatomical landmarks recorded on microCT scans of Apert syndrome mice at P0. (B) Autopod (hand). Landmarks were recorded at the midpoint of the proximal and distal tips of the distal, mid and proximal phalanges (P1–P28) and the metacarpals (P29–P38). Proximal phalanx I, middle phalanx V, and metacarpal I are not displayed because these bones have not yet mineralized at P0 and could not be visualized in many specimens. (C) Zeugopod. Landmarks were recorded at the midpoint of the proximal and distal tips of the radius (R1–R2) and ulna (U1–U2). (D) Stylopod. Landmarks were recorded at the midpoint of the proximal and distal tips of the humerus (H1–H2), as well as at the tip of the deltoid process (H3) (Figure 1—source data 1). (E) Scapula. Landmarks were recorded at the most superior and inferior lateral points of the scapula (S1–S2), the most posterior point of the spine (S3), the most antero-medial point of the acromion process (S4) and the medial, superior and inferior points of the glenoid cavity (S5–S7) (Figure 1—source data 2). (F) Clavicle. Landmarks were recorded at the medial point of the sternal and the acromial ends (C1–C2). (G) Length and volume differences in the forelimbs of Apert syndrome mouse models. Schematic representation of the forelimb of a P0 mouse showing, in different colors, statistically significant differences in bone length and volume (as measured by two-tailed one-way ANOVA or Mann-Whitney U-test) in Fgfr2+/P253R mice and unaffected littermates, as specified in Table 1. Longer/shorter refers to length, whereas larger/smaller refers to volume. (H, I) Shape differences in the forelimbs of Apert syndrome mouse models. Scatterplots of PC1 and PC2 scores based on Procrustes analysis of anatomical landmark locations representing the shape of the left humerus (H) and the left scapula (I) of unaffected (N = 10) and mutant (N = 12) littermates of Apert syndrome mouse models (Figure 1—source data 1 and Figure 1—source data 2). Convex hulls represent the range of variation within each group of mice.

Figure 2.. Qualitative visualization of Dusp6 gene expression in unaffected and Fgfr2+/P253R mouse embryos at E10.5 and E11.5.

OPT scans of embryos WMISH stained for Dusp6 revealed the anatomical location of Dusp6 gene expression (shown in yellow). For each stage, the main expression domains are highlighted on the left for anatomical reference on a lateral view of a 3D reconstruction of a Fgfr2+/+ unaffected embryo (fb, forebrain; fl, forelimb; hb, hindbrain; hl, hindlimb; ht, heart; is, isthmus; lnp, lateral nasal process; md, mandibular prominence; mnp, medial nasal process; mx, maxillary prominence; so, somites). On the right, 3D reconstructions of three unaffected and three Fgfr2+/P253R mutant embryos from the same litter are displayed to represent the high degree of variation in developmental age within litters. Embryos are not to scale. Original 2D images of the Dusp6 WMISH experiments are available in Figure 2—source data 1.

Figure 3.. New quantitative analysis method for 3D gene expression data, based on Geometric Morphometrics.

Mouse embryos between E10.5 and E11.5 were analyzed with WMISH to reveal the expression of Dusp6 (1), and then cleared with BABB, and OPT scanned using both fluorescence and transmission light (2). The external surface of the embryo was obtained from the 3D reconstruction of the fluorescence scan (2). Multiple thresholding of the transmission scan by choosing different levels of grey values as shown by the histogram allowed the visualization of gene expression patterns at different intensities (3). Moderate gene expression (shown in green) was displayed as the isosurface obtained using a threshold of the grey value computed as 2/3 of the last grey value showing the Dusp6 expression domain (3). High gene expression (shown in yellow) was displayed as the isosurface obtained using a threshold of the grey value computed as 1/3 of the last grey value showing the Dusp6 expression domain (3). From the whole mouse embryo isosurfaces, all four limb buds were segmented (5A). From the high gene expression isosurface, the Dusp6 domains from all of the available limbs were segmented (5B). Maximum curvature patterns were displayed to optimize landmark recording (5). For each limb, we captured the shape and size of the limb bud (6A) and the underlying high Dusp6 gene expression pattern (6B), recording the 3D coordinates of anatomical landmarks (yellow dots), curve semi-landmarks (red dots) and surface semi-landmarks (pink dots). Anatomical and curve landmarks were recorded manually on each limb. Surface landmarks were recorded on one template limb and interpolated onto target limbs (Video 1). Landmark coordinates were the input for Geometric Morphometric (GM) quantitative shape analysis (7, 8) that was used to superimpose the landmark data (GPA, Generalized Procrustes Analysis), to compute limb size (centroid size), and to explore both shape variation within limbs and gene expression domains within litters by Principal Component Analysis (PCA). Finally, the covariation patterns for the shape of the limb and the shape of the gene expression domain were also explored using the Partial Least Squares (PLS) method.

Figure 4.. Developmental variation within litters of the Apert syndrome mouse model.

All of the limbs were individually staged using our publicly available web-based staging system (https://limbstaging.embl.es). The stage of the limb bud was estimated after alignment and shape comparison of the spline with an existing dataset with a reproducibility of ±2 hr. As both Litter 1 and Litter 2 comprised limbs between E10 and E10.5, we pooled these limbs into the 'Early' period. Litter 3, comprising limbs between E10 and E11, represented the 'Mid' developmental period. And finally, Litter 4, which included limbs between E11 and E11.5, was considered as the 'Late' period. In subsequent analyses, forelimbs and hindlimbs were analyzed separately, as specified in Table 2. Increasing color intensities in dots represent an increasing number of limbs with the same staging result. Individual scores are available in Figure 4—source data 1.

Figure 5.. Tracing of limb phenotypes (anatomical and molecular) back through developmental time to the earliest moment of appearance.

Principal Component Analysis based on the Procrustes-based semi-landmark analysis was used to assess the shape of the limbs and the corresponding Dusp6 expression domains for each developmental period (Figure 5—figure supplement 1). Each period was analyzed separately for the shape of the limb (A–D) and the Dusp6 expression domain (E–H), as specified in the 'Materials and methods' and in Table 2. Scatterplots of PC1 and PC2 axes with the corresponding percentages of total morphological variation explained are displayed for each analysis, along with the morphings associated with the negative, mid and positive values of the PC axis that separates mutant and unaffected littermates (PC1 or PC2, as highlighted in bold black letters in the corresponding axis). Morphings are displayed in grey tones when the analysis showed no differentiation between mutant and unaffected littermates. Morphings are displayed in color when the analysis revealed differentiation between mouse groups (blue: unaffected littermates; pink: mutant littermates). Limb buds and Dusp6 domains are not to scale and are oriented with the distal aspect to the left, the proximal aspect to the right, the anterior aspect at the top and the posterior aspect at the bottom of all images. Convex hulls represent the ranges of variation within each group of mice. Individual scores are available in Figure 5—source data 1. Results from the complete dataset analyzing hindlimbs (Figure 5—figure supplement 2, Figure 5—figure supplement 2—source data 1) and forelimbs (Figure 5—figure supplement 3, Figure 5—figure supplement 3—source data 1) from each litter separately are also available.

Figure 5—figure supplement 1.. Anatomical landmarks recorded on the 3D reconstructions of the limbs and the Dusp6 expression domains obtained from OPT scans of Apert syndrome embryos at E10.5–E11.5.

In the limb, we manually recorded the 3D coordinates of five anatomical landmarks (A): L1–L2, the most medial points on the anterior and the posterior sides of wrist, collected along the curve that divides the limb on the dorsal and ventral sides; L3, the most distal point of the limb, collected along the curve that divides the limb on the dorsal and ventral sides; and L4–L5, the most superior and inferior points on the dorsal and ventral sides of the limb. We also recorded between 25 and 60 points along the curve that divides the limb on the dorsal and ventral sides to define an outline spline of the limb. Using Viewbox 4, we obtained the 3D coordinates of 19 equidistant points along the outline spline that were considered in further comparative shape analyses to be curve semi-landmarks (LC1–LC19) (B). Finally, we recorded 52 surface semi-landmarks on the dorsal and ventral sides of the limb interpolating the points from a reference limb to each target limb, using the five anatomical landmarks manually recorded as matching points and by minimizing the bending energy (LS1–LS52) (C). In the gene, we manually recorded the 3D coordinates of four anatomical landmarks (D): G1–G2, the points at the most anterior and posterior tips of the gene expression domain; G3-G4, the points at the maximum curvature of the proximal side of the gene expression domain, on the dorsal and the ventral sides (D). (We also recorded 10–30 points along the distal and proximal curves of the gene expression domain at the sagittal plane, and 30–60 points along the outline of the gene expression domain, on its proximal side.) Using Viewbox 4, we obtained the 3D coordinates of eleven equidistant points along the spline of the distal and proximal curves of the gene expression domain (GCd1–GCd11 and GCp1–GCp11) (E and F), and 19 equidistant points along outline of the gene expression domain, on its proximal side (GCo1–GCo19) (G) that were considered in further comparative shape analyses as curve semi-landmarks, and are displayed simultaneously in (H).

Figure 5—figure supplement 2.. Principal Component analyses based on the Procrustes-based semi-landmark analysis of the shape of the hindlimbs and the corresponding Dusp6 expression domains for each developmental period.

Scatterplots of PC1 and PC2 axes with the corresponding percentage of total morphological variation explained are displayed for each analysis, along with morphings associated with the negative, mid and positive values of the PC axis that separates mutant and unaffected littermates (PC1 or PC2, as highlighted in bold black letters in the corresponding axis). Morphings are displayed in grey tones when the analysis shows no differentiation between unaffected and mutant littermates. Morphings are displayed in color when the analysis reveals differentiation between the mouse groups (blue: unaffected littermates; pink: mutant littermates). Limb buds and Dusp6 domains are not to scale and are oriented distally to the left, proximally to the right, anteriorly to the top and posteriorly to the bottom. Convex hulls represent the ranges of variation within each group of mice, as defined in Table 2 and Figure 4. Individual scores are available in Figure 5—figure supplement 2—source data 1.

Figure 5—figure supplement 3.. Principal Component analyses based on the Procrustes-based semi-landmark analysis of the shape of the forelimbs and the corresponding Dusp6 expression domains for each developmental period.

Scatterplots of PC1 and PC2 axes with the corresponding percentage of total morphological variation explained are displayed for each analysis, along with morphings associated with the negative, mid and positive values of the PC axis that separates mutant and unaffected littermates (PC1 or PC2, as highlighted in bold black letters in the corresponding axis). Morphings are displayed in grey tones when the analysis shows no differentiation between unaffected and mutant littermates. Morphings are displayed in color when the analysis reveals differentiation between the mouse groups (blue: unaffected littermates; pink: mutant littermates). Limb buds and Dusp6 domains are not to scale and are oriented distally to the left, proximally to the right, anteriorly to the top and posteriorly to the bottom. Convex hulls represent the ranges of variation within each group of mice, as defined in Table 2 and Figure 4. Individual scores are available in Figure 5—figure supplement 3—source data 1.

Figure 6.. Quantitative correlation between the size and shape of the limbs and the Dusp6 expression pattern.

(A–D) Comparison of limb bud size and Dusp6 volume in unaffected and Fgfr2+/P253R mutant littermates across development (Table 3). Limb size was measured as limb centroid size (A), whereas the size of the Dusp6 expression was measured as the volume of the gene domain (B), as specified in the 'Materials and methods' and in Table 2. Statistically significant differences as revealed by two-tailed t-tests are marked with asterisks, representing the degree of significance: *P-value=0.03, **P-value=0.01. Results from hindlimbs and forelimbs are separately available in Figure 6—figure supplement 1. The association between the size of the limbs and the volume of the Dusp6 domain was assessed separately for forelimbs (C) and hindlimbs (D). Individual scores for all of these analyses are available in Figure 6—source data 1. (E) Time course assessing the morphological integration pattern between the limb phenotype and the shape of the gene expression pattern using partial least squares analyses. Associated shape changes from late to early limb development are shown from morphings associated with the negative, mid and positive values of PLS1, which accounted for almost all of the covariation (97.6% in forelimbs and 99.5% in hindlimbs) between the limb buds and the Dusp6 gene expression domains (Figure 6—figure supplement 2). Limb buds and Dusp6 domains are not shown to scale and are oriented distally to the left, proximally to the right, anteriorly to the top and posteriorly to the bottom. On the right, representing the positive extreme of PLS1 axis, typical early limb buds showed a protruding shape (i.e. short in the proximo-distal axis and symmetrical in the antero-posterior axis) associated with a flat-bean shaped Dusp6 expression zone localized in the distal limb region, underlying the apical ectodermal ridge and spreading proximally towards the dorsal and ventral sides of the limb. On the left, representing the negative extreme of the PLS1 axis, limb buds were elongated in the proximal axis and asymmetric on the antero-posterior axis, with an expansion of the distal limb region and a contraction of the proximal region, at the wrist level. This limb shape, which is typical of more developed limbs, was associated with a Dusp6 expression that was extended underneath the apical ectodermal ridge towards the anterior and the posterior ends of the gene expression zone, but reduced on the dorsal and ventral sides of the limb. Individual scores for all of these analyses are available in Figure 6—figure supplement 2—source data 1. Legends for supplementary figures.

Figure 6—figure supplement 1.. Quantitative size analyses of limb bud and Dusp6 volume in unaffected and Fgfr2+/P253R mutant littermates throughout development.

Comparison in hindlimbs (A, C) and forelimbs (B, D) at each developmental period, as defined in Table 2 and Figure 4. Limb size was measured in terms of limb centroid size (μm), whereas the size of the Dusp6 expression zone was measured as the volume of the gene expression domain (μm3). Statistical significant differences as revealed by two-tailed t-tests are marked with asterisks: *P-value<0.05. Individual scores are available in Figure 6—source data 1.

Figure 6—figure supplement 2.. Morphological integration between the shapes of the limb bud and the Dusp6 gene expression domain, as measured by partial least square (PLS) analysis.

Results display the morphings associated with the negative, mid and positive values of PLS1 for hindlimbs (A) and forelimbs (B), which accounted for more than 95% of the covariation between the limb buds and the Dusp6 gene expression domains. Separate PLS analyses for hindlimbs (C) and forelimbs (D) for each developmental group, as defined in Table 2 and Figure 4, are also shown. Formal testing of differences in the intensity of the integration patterns in unaffected and Fgfr2+/P253R mutant mice at each stage was not possible because of sample size limitations, as specified in 'Materials and methods' and in Table 2. Individual scores are available in Figure 6—figure supplement 2—source data 1. Rich media.