The Drosophila blastoderm embryo undergoes rapid cycles of nuclear division. This poses a challenge to genes that need to reliably sense the concentrations of morphogen molecules to form desired expression patterns. Here we investigate whether the transcriptional state of hunchback (hb), a target gene directly activated by the morphogenetic protein Bicoid (Bcd), exhibits properties indicative of inheritance between mitotic cycles. To achieve this, we build a dataset of hb transcriptional states at the resolution of individual nuclei in embryos at early cycle 14. We perform a spatial point pattern (SPP) analysis to evaluate the spatial relationships among the nuclei that have distinct numbers of hb gene copies undergoing active transcription in snapshots of embryos. Our statistical tests and simulation studies reveal properties of dispersed clustering for nuclei with both or neither copies of hb undergoing active transcription. Modeling of nuclear lineages from cycle 11 to cycle 14 suggests that these two types of nuclei can achieve spatial clustering when, and only when, the transcriptional states are allowed to propagate between mitotic cycles. Our results are consistent with the possibility where the positional information encoded by the Bcd morphogen gradient may not need to be decoded de novo at all mitotic cycles in the Drosophila blastoderm embryo.
Citation: He F, Ma J (2013) A Spatial Point Pattern Analysis in Drosophila Blastoderm Embryos Evaluating the Potential Inheritance of Transcriptional States. PLoS ONE 8(4): e60876. doi:10.1371/journal.pone.0060876
Editor: Amy S. Gladfelter, Dartmouth College, United States of America
Received: January 8, 2013; Accepted: March 4, 2013; Published: April 9, 2013
Copyright: © 2013 He, Ma. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Funding: This work was supported in part by a fellowship award from American Heart Association (10POST4170015) to FH and by grants from National Institutes of Health (1R01GM101373) and National Science Foundation (IOS-0843424) to JM. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors have declared that no competing interests exist.
Transcription is a molecular process where the enzyme RNA polymerase (RNAP) makes RNA copies from the DNA templates –. It is an inherently stochastic process that is reflective of the random nature of the underlying molecular events necessary for RNAPs to transcribe through a gene. These molecular events may include, for example, the remodeling of chromatin, the binding or unbinding of transcription factors and the loading or clearance of RNAPs –. The prevailing theoretical framework on regulation of transcription considers only the “recent” biochemical events that lead to successful transcription of a gene –. But can a gene’s transcriptional state be also influenced by the “past” events or its “past” experiences? We raise this question because mitosis is known to abort transcription , , yet active transcription of many patterning genes becomes detectable almost immediately after a mitotic cycle in the Drosophila blastoderm embryo –. Thus the discussions about the recent and past events in our current work make an exclusive reference to mitosis.
In Drosophila, embryonic patterning along the anterior-posterior (AP) axis requires the morphogen gradient of Bicoid (Bcd) –. Bcd is a transcriptional activator that binds to the enhancers of its target genes such as hunchback (hb) and activates their transcription , , , . The mechanisms of sensing the concentrations of morphogen molecules by a gene or cell are crucial to our knowledge of how the positional information provided by a morphogen gradient is decoded during development –. Theoretical studies suggest that there are fundamental limits to the time period necessary for a cell or gene to accurately “sense” or “read” the differences in morphogen concentrations to make a response to form a desired expression boundary , . This time period is dependent on, among other things, both the diffusion constant and the nuclear concentration of Bcd at the boundary position. Based on available measurements, the time it takes the hb gene to accurately read the nuclear Bcd concentration and form a de novo expression boundary has been calculated to be on the order of tens of minutes or longer , , . This is in contrast to the observation that active transcription (at and beyond the boundary) of hb gene copies can take place almost immediately after mitosis, on the order of less than a minute , . To reconcile these conflicting properties of the early Drosophila embryos, a memorization mechanism of an unknown nature has been proposed, where the hb transcriptional state of a nucleus can be influenced by the state of its parent nucleus . This hypothesis provides a rational explanation for how the positional information provided by the Bcd gradient can be decoded to form reliable and precise boundaries of its target genes in early embryos that are undergoing rapid mitotic cycles in roughly every 10 minutes, a time period that would not allow de novo decoding and boundary formation. However, due to a lack of experimental tools suitable for tracking the transcriptional states of hb gene copies between mitotic cycles, this important hypothesis has remained untested.
In this study, we perform a spatial point pattern (SPP) analysis to evaluate the spatial properties of hb transcriptional states at the resolution of individual nuclei. SPP relates to the distribution characteristics of, i.e., spatial relationships among, a series of point locations . Many biological examples, such as a sheet of cells, can be reduced to a pattern of mapped points , . Here we treat each nucleus within the monolayer of a blastoderm embryo as a point and its transcriptional state as an event. Our study suggests properties of spatial clustering for nuclei with both or neither copies that are undergoing active transcription in snapshots of embryos. Our statistical and simulation analyses are consistent with the possibility that the hb transcriptional states of the nuclei can be propagated (i.e., inherited) between mitotic cycles. The statistical and simulation methods used in our current study may be of general value in investigating developmental problems that require knowledge about the spatial properties of cellular and molecular decisions.
Materials and Methods
Our dataset was constructed based on our FISH data detecting the nascent transcripts of hb in early cycle 14 wildtype embryos (w1118) at the peak time of active transcription of hb gene copies , . These nascent transcripts were detected by a 283-nucleotide probe against the hb intron RNA sequence. This intron is located 145 bp downstream of the transcription start site from the Bcd-responsive hb promoter P2. Thus, the detected signals, which were discrete fluorescent dots inside the nucleus, referred to as intron dots, capture individual copies of hb undergoing active transcription near the P2 promoter , , . Our current dataset consists of data from 14 whole-mount embryos, which had been flattened to maximize the number of nuclei captured. In our analysis, ~6 confocal z-section images at 0.5 µm intervals were taken to capture all intron dots in the nuclei in one of the flattened sheets of the nuclear monolayer. Thresholds to identify an intron dot included both signal cluster size and signal intensity, and they were optimized to minimize false positive and false negative cases and to detect a stable intron dot expression pattern along AP axis as discussed previously .
Each embryo was oriented with its major axis (i.e., AP axis) parallel to the x-axis on the image, with the anterior tip set as the origin. The center of each nucleus was identified by MATLAB and its x-y coordinates were recorded. Each image contains 3000~4000 nuclei within the area of ~600 µm×250 µm (Figures 1 and S1). An experimental field of 100 µm×100 µm (i.e., 130 µm to 230 µm on x-axis and −50 µm to 50 µm on y-axis) was cropped for further analysis (Figures 1 and S1). We chose this location after considering the following specific issues. 1) We avoided the edges of the embryos to minimize errors due to out-of-focus detection of nuclei. 2) We avoided the anterior tip of the embryo, where the nuclear density is known to deviate significantly from other parts of the embryo . 3) We avoided the boundary of the hb expression domain because this boundary itself represents a spatial pattern that would interfere with our analysis. Changing the location of the cropped field by ±5 µm along either axis did not affect any of our conclusions.
Shown is an embryo with nuclei masked according to hb transcriptional states, i.e., the number of hb intron dots detected. Blue: 0-dot nuclei; green: 1-dot nuclei; red: 2-dot nuclei; black: nuclei with over 2 dots. (B) is the cropped experimental field from the solid square in (A). See Figure S1 for locations of the experimental fields in all 14 embryos.
Characterizations of Nuclear Density and Inter-nuclear Distances
The overall densities of the nuclei of individual types or all types combined were calculated as d = N/A, where N denotes the number of nuclei of the type of interest within the entire field that has an area size A = 100×100 µm2 (see Results for definitions of nuclear types). The local densities of the nuclei (of either individual or combined types) were estimated at each nuclear position in the pattern by di = ni/πR2, where ni denotes the number of nuclei of a type within a distance R from the center of nucleus i. Here we chose R = 8 µm, a value greater than the diameter of a typical nucleus (~6 µm) at this stage of embryonic development.
Center-to-center distances between all pair-wise nuclei were computed as follows. A nucleus in an experimental field was identified as the reference and the distance from that nucleus to every other nucleus in the field was measured. This procedure was repeated for every nucleus as a reference on this field and then repeated for each of the 14 experimental fields. For a given nucleus, the nearest-neighbor distance (NND) is the minimum value among its distances to all other nuclei . We calculated the regularity index (RI), which was defined as the mean NND divided by its standard deviation . Theoretically, for a random distribution of points in two dimensions, RI ≈ 1.91 . This value was used as a standard against which the RI values obtained from our experimental or simulated fields were evaluated (Figure 2).
Shown are regularity indexes (RI) derived by the measurement of nearest-neighbor distances. Black: results from 14 experimental fields; white: results from 14 random simulations, each of which was constrained by the total nuclear number and the proportions nuclear types observed in the corresponding experimental field. The dashed line shows the theoretical expectation of 1.91 for randomly distributed points (see also Materials and Methods). In simulations, nuclei of different types or all nuclei combined behave indistinguishable from each other or the theoretical expectation. Error bars show the standard deviations calculated from 14 experimental or simulated fields. p-values from paired Student’s t-tests to compare the means of experimental and simulated RI values were shown over the bars.
Since our experimental embryos were flattened to maximize the number of intron dots captured, we performed an analysis to assess the impact of embryo flattening on the spatial distributions of nuclei. For each blastoderm embryo (N = 10), we performed confocal imaging both before (i.e., free in mounting medium) and after flattening. We found that flattening an embryo led to an expansion of ~30 µm along x- or y-axis (from 562±17 µm to 590±20 µm along x-axis and from 224±21 µm to 255±22 µm along y-axis). The distances between neighboring nuclei in the flattened sheets are increased by ~4% relative to the 3-dimentional distances between the neighboring nuclei on the periphery of non-flattened embryos. Under both conditions, the nuclei were found to be regularly and evenly spaced, suggesting that flattening did not cause any significant localized nuclear movements within the experimental fields chosen for our analysis.
To calculate Ripley’s K-functions , we binned the measured distances and counted the overall number of nuclei at each distance r from a reference nucleus. In effect, our exercise created annuli around the reference nucleus and measured the number of nuclei inside circles at different r intervals. “Edge effect” may arise for a reference nucleus within distance r from the boundaries of a field. We compensated this effect by surrounding the field under examination with 8 identical copies of this field (effectively forming a larger sheet that contains 9 identical “tiles”). By doing this, we were able to analyze the complete SPPs of any nuclei within the field (the center tile), including those close to the edge. We used the 2 µm r intervals for the annuli in our current study; we also performed analyses with increasing interval sizes and observed comparable results up to a 5 µm r interval size. Eq. 1 defines K-function, K(r), where d is the nuclear density.(1)
In a bivariate pattern, the nuclear densities of Type i and Type j are di and dj, respectively. A complete description of the second-order properties of the bivariate pattern requires a consideration of all possible pairs between the opposite types. We define the bivariate K-functions K ij(r) by(2)
Various benchmark hypotheses when described by the K-function have well-defined expressions. In particular, the CSR hypothesis states that K(r) = πr2. Independence between two opposite types has Kij(r) = πr2, irrespective of the individual point patterns or the forms of Ki(r) and Kj(r). In our study, all K-functions were transformed to H-functions to allow visual illustration of departures from the benchmark hypotheses (see Eqs. 3–4 and the related texts in Results).
Tests of Significance
To test statistical significance, we used a Monte Carlo approach for estimating two properties of each experimental field, 1) the K-functions and 2) the proportions of different types of nuclei that surround a nucleus of a given type. In this simulation, we used each experimental field as a “mold” for the locations of all the nuclei inside the field. Then, we randomly assigned a type to each and every nucleus while maintaining the overall proportions of the three types as observed experimentally for that mold. To evaluate significance, we performed Kolmogorov-Smirnov tests to compare experimental and simulated data.
To consider the proportions of different types of nuclei that surround a nucleus of a given type, we used the Voronoi diagram that identifies the immediate-neighboring nuclei that share an edge in a mold . For each and every nucleus, we recorded the number of immediate-neighboring nuclei belonging to each of the three types in either experimental or simulated fields. For each mold, we performed 199 simulations and estimated the 95% confidence interval (CI) using the results that ranked as the top or bottom 5 simulations. In addition to the Monte Carlo approach, we also theoretically predicted the CI by treating each of the 14 experimental patterns as a draw from a probability distribution. Specifically, we assumed that, in a population of N nuclei, the number of a given type (n) follows a binomial distribution with probability Pi, i.e., n ~ B(N, Pi). Then the 95% CI was determined by(3)The parameter Pi was estimated for each nuclear type in the experimental fields.
Modeling Nuclear Lineages and Mitotic Inheritance
To simulate the nuclear lineages, we used the nuclear positions shown in Figure 1B as the mold and randomly assigned every 8 neighboring nuclei within the mold into a clone. This clone represents the nuclei that would have been derived from a single nucleus at cycle 11 after undergoing 3 mitotic cycles. Then, we randomly divided each of the 8-nuclei clones into two 4-nuclei clones. These two clones represent the nuclei that would have been derived from the two nuclei at cycle 12 after undergoing 2 mitotic cycles. Finally, we randomly divided each of the 4-nuclei clones into two 2-nuclei clones, which represent the nuclei that would have been derived from two nuclei at cycle 13 after undergoing one mitotic cycle. Effectively, this exercise led to the generation of nuclear lineages from cycle 11 to 14 without specifically modeling or knowledge of the nuclear positions and movements prior to cycle 14.
To simulate the spatial patterns of transcriptional states of nuclei in blastoderm embryos, we assumed, based on reported experimental data , that the probability of hb gene copies undergoing active transcription (referred to as the transcribing probability thereafter) at cycles 11 to 14 is 90%, 90%, 50% and 50%, respectively. In simulations where there is no inheritance, we assumed the number of active copies in each nucleus as a random draw from a binomial distribution, which is a function of transcribing probability at each cycle. In simulations where inheritance was allowed, the two daughter nuclei (at each new cycle) were assigned to exist as the same type as their parental nucleus. Upon entering the cycle 13 interphase, active copies were randomly switched off to achieve the experimentally observed reduction in the transcribing probability. We performed 99 simulations for each of the two scenarios.
The SPP Dataset and Basic Features of the Nuclear Monolayer
To analyze the spatial point pattern (SPP) of nuclei with distinct hb transcriptional states, we constructed a dataset from images of 14 Drosophila embryos at early nuclear cycle 14 (see Materials and Methods). These images represent snapshots capturing robust hb transcription at the resolution of individual copies of the gene inside the nuclei of the embryo . Here an intronic probe for hb was used in FISH to detect its nascent transcripts in whole-mount embryos. These nascent transcripts were detectable as distinct fluorescent dots, referred to as intron dots, each representing, in the snapshot, a copy of hb undergoing active transcription near the Bcd-responsive P2 promoter (see Materials and Methods). Our dataset includes the identities and locations of the nuclei containing 0, 1 or 2 intron dots (Figure 1A blue, green and red). We define all nuclei with a same number of intron dots as being one type. Thus there are three types of nuclei: Type-0, Type-1 and Type-2 (a minimal number of nuclei with erroneously more than 2 dots detected were also treated as Type-2 nuclei; Figure 1A black). The fractions of the three types of nuclei in our current dataset (see Table S1 for nuclear counts in individual embryos) are well explained by a binomial distribution, which argues against bistability of transcriptional states .
It has been shown that Bcd induces hb transcription in a highly cooperative manner . This cooperative action ,  effectively divides an embryo into two broad domains along the anterior-posterior (AP) axis , , . Here we refer to the anterior portion of the embryo with Bcd-dependent active hb transcription as the transcribing domain, whereas the posterior portion lacking Bcd-dependent active hb transcription as the non-transcribing domain. In our current study, we focus exclusively on the transcribing domain. Here our goal is not to analyze directly how the two domains are established (i.e., boundary formation as a function of the Bcd input in a threshold response) but, rather, within the transcribing domain, how the nuclei with distinct hb transcriptional states are distributed spatially. To facilitate our analysis, we cropped the images to obtain areas (i.e., experimental fields) of an identical size within the transcribing domains of individual embryos (see also Materials and Methods and Figure 1B). Since visual inspections (Figure 1B) are not adequate for extracting and describing the characteristics of the spatial relationships among the nuclei of a given type, we approached this problem through the use of statistical methods.
Several methods have been developed to analyze SPP data , . In our study, we used first-order statistics to characterize the basic features of the nuclear monolayer in cycle-14 embryos with regard to nuclear density and variations. Here the experimental fields have a mean nuclear density (d) of 0.024 nuclei/µm2 (averaged among different embryos) and a standard deviation (SD) of 0.002 nuclei/µm2. The overall densities of Types-0, -1 and -2 nuclei in the fields are 0.004±0.001, 0.011±0.002 and 0.009±0.002 nuclei/µm2, respectively. Previous studies have shown that an embryo at early cycle 14 interphase has relatively uniform nuclear density throughout, with a dominant exception at its terminal regions , , . In our chosen experimental fields (see Materials and Methods), the local densities of all nuclei, Types-0, -1 and -2 nuclei are 0.026±0.003, 0.013±0.005, 0.013±0.004 and 0.015±0.004 nuclei/µm2, respectively (all values are averaged among different nuclear positions with SD shown). These results show that, consistent with previous reports, our experimental fields, which were chosen to be devoid of the terminal regions of the embryo (see Materials and Methods and Figure S1), do not exhibit major variations in local nuclear densities.
To further evaluate the basic features of the nuclear monolayer, we calculated the regularity index (RI). RI is defined as the ratio of the mean of the nearest-neighbor distance (NND) to its standard deviation. It measures the uniformity of NND, providing a widely-used criterion for evaluating spatial regularity ,  (see also Materials and Methods). NND calculation does not take into account the actual size of the nucleus. Thus the relatively uniform size and spacing of the nuclei in the monolayer are features of spatial regularity of the nuclei themselves under the framework of the RI analysis (see also below for further support to this notion). Figure 2 plots the RI values of the three types of nuclei and all types combined (solid bars; RI = 2.81±1.61, 4.49±1.29, 3.90±1.48 and 12.49±2.96, for Types-0, 1, 2, and all types, respectively). To statistically evaluate spatial regularity, we generated a random distribution of nuclei within a two-dimensional field and calculated RI values for the simulated data. Here, both the field size and nuclear number, but not nuclear positions, were constrained by our experimental measurements of an experimental field. This analysis was performed for each of the 14 experimental fields. We found that the experimentally measured RI values (Figure 2 solid bars) are significantly higher than their corresponding values derived from random simulations (Figure 2 white bars; RI = 1.89±0.14, 2.01±0.17, 1.87±0.09 and 1.97±0.11, respectively; p = 0.04, 10−7, 10−5 and 10−13, respectively). In each case, these values are also higher than the theoretical value calculated for random spatial distributions of points (see Materials and Methods for details about this theoretical value and simulations). These results confirm quantitatively that the nuclei themselves within our experimental fields of the nuclear monolayer have regular spacing and uniform size.
Type-0 and Type-2 Nuclear Patterns are Formed Non-randomly
To investigate the spatial distributions of the nuclei, we conducted an analysis using second-order statistics. Unlike the distance-based NND method, second-order statistics describes the characteristics of the point distributions as a function of area size. Here we used Ripley’s K-function, which has been extensively employed in analyzing spatial patterns in biology , . The K-function analysis uses all inter-nuclei distances to evaluate a spatial point pattern. It measures the number of events that occur inside the area of a circle with a radius r from a reference event. In essence, a circle of a specific radius r is passed over a pattern, and the number of events that occur inside this circle is counted. This counting is repeated at different radii for this nucleus and then repeated for all nuclei. Thus, the K-function is presented as a plot against r. Note that the definition of K-function has been arranged such that it is only a function of r. It is independent of the nuclear density d of the experimental fields (see Materials and Methods). Since the K-function analysis is conducted at multiple radii, it can reveal a pattern’s characteristics that are dependent on the scales of neighborhood size being evaluated. This particular feature of K-function analysis circumvents the limitations of the NND method and, thus, can reveal deeper information about spatial patterns.
Subtraction of r facilitates a visual presentation illustrating departures from complete spatial randomness (CSR). In this analysis, CSR is the null hypothesis stating that the positioning of the nuclei of a given type (the events) within the experimental fields follows a homogeneous Poisson process. Thus violations of this null hypothesis may arise from major variations in the local density of events. In Eq. 4, a positive Hi(r) value means that, on average, there are more nuclei of Type i inside the area of a circle with a radius r than expected of CSR. Thus, a positive Hi(r) value indicates spatial clustering when evaluated at r. Conversely, a negative Hi(r) value indicates repulsion among the nuclei at r.
Figure 3 shows the H-functions of different types of nuclei plotted against r. On scales approaching the nuclear diameter (4~8 µm at early cycle 14), the H-functions of the nuclei of each individual type or all types combined are all negative, indicative of a repulsive effect among these nuclei. This result is simply reflective of the basic feature of the nuclear monolayer discussed above, namely, the relatively uniform size and spacing of the nuclei themselves. However, on scales exceeding the typical size of the nucleus, the H-functions of Types-0 and -2 nuclei both become positive (Figure 3A and C). These results show that the nuclei of these two types have properties of spatial clustering. In our analysis, Type-1 nuclei exhibit no significant difference from CSR (Figure 3B).
Shown are the profiles of H-functions with respect of the distance r for Type-0 nuclei (A), Type-1 nuclei (B), Type-2 nuclei (C), Type-0/Type-2 bivariate (purple dots in D) and all nuclei combined (solid lines in D). The dashed line shows the theoretical expectation of complete spatial randomness (CSR; see also Materials and Methods). Error bars show the standard deviations calculated from 14 experimental fields.
Type-0 and Type-2 Nuclei are Distributed as Dispersed Clusters
If Types-0 and -2 nuclei are indeed spatially clustered, there is a possibility that the nuclei of these two opposite types may exhibit properties of mutual repulsion. To evaluate this possibility, we calculated the bivariate H-function (Eq. 5). Here we used a nucleus of one type as a reference point (i.e., Type i) and counted the number of nuclei of the opposite type (i.e., Type j, i ≠ j) that fall within the area of a circle with radius r.(5)
For bivariate point patterns, where patterns of two opposite events (i.e., the nuclei of two different types) are considered in relation to each other, the null hypothesis is independence. Under this null hypothesis, Hij(r) = 0. Hij(r) >0 indicates that, on average, there are more Type-j nuclei (than expected of independence) inside the area of a circle with a radius r from a Type-i nucleus. This suggests an attraction between the nuclei of the opposite types. Conversely, Hij(r) <0 suggests repulsion between the nuclei of the opposite types. Either attraction or repulsion would reject the null hypothesis stating that these two opposite events are determined independently of each other. Based on the bivariate H-function analysis (Figure 3D), Types-0 and -2 nuclei exhibit properties of mutual repulsion (Figure 3D purple). Using the same approach, the null hypothesis cannot be rejected when the two opposite events are Type-0 and Type-1 nuclei or are Type-1 and Type-2 nuclei (data not shown).
To extend the evaluation and understanding of our experimental data, we generated randomly simulated data based on each of the “molds” derived from the 14 experimental fields (see Materials and Methods). Figure 4 plots the H-function profiles of the simulated data. Here we used Kolmogorov-Smirnov tests to further evaluate the differences between the H-functions of the experimental and simulated data (shown in Figures 3 and 4, respectively). For each pattern, we determined the range of r, within which the H(r) differences are significant (p<0.05). The r ranges for Type-0 analysis, Type-2 analysis and Type-0/Type-2 bivariate analysis are 6~44 µm, 12~80 µm and 6~100 µm. They document that the experimental and simulated data differ in a consistent manner within broad r ranges, which contrasts with the narrow r range of 30~34 µm for Type-1 analysis. Together these results provide further support to the observation that, in experimental fields, Types-0 and -2 nuclei are spatially clustered and these two types of nuclei exhibit mutual repulsion.
Shown are results of an analysis that is identical to that of Figure 3, but simulated fields are used here. Unlike simulations shown in Figure 2, the simulations performed here were based on the molds derived from the experimental fields, with both the locations of the nuclei and the proportions of the nuclear types constrained by experimental observations. Note that the profiles of individual nuclear types are all below 0. This is expected because the nuclei of a given type are spatially mixed with the nuclei of other types when randomly distributed and thus the nuclei of a given type always exhibit a self-repelling SPP. Moreover, because the proportion of Type-0 nuclei is the smallest, the H-function profile of this type is the most negative among all three types. Also note that the profiles do not converge to zero at r = 100 µm, which is the boundary of the reference nuclei to compute distance, but not the limit of a complete SPP (see also Materials and Methods).
Distributions of Different Types of Nuclei in Spatial Relations to Type-2 Nuclei
If Type-2 nuclei in our experimental fields form local clusters and, conversely, Type-0 and Type-2 nuclei exhibit mutual repulsion, there is a corollary that can be tested. It states that, from a given Type-2 nucleus, it would take a shorter distance to reach another Type-2 nucleus than to reach a Type-0 nucleus in our experimental fields. To test this corollary through the use of an independent method (other than the use of Ripley’s K-function), we plot the empirical cumulative distribution functions for the distances between a Type-2 nucleus and another nucleus of any type. Figure 5A shows that pairs of Type-2 nuclei are more likely to be closer (blue) than pairs of nuclei of different types (green and red). Kolmogorov-Smirnov tests suggest that, at the 1% significance level, all these three distributions are different.
Shown are the empirical cumulative distribution function plots of the distances from a Type-2 nucleus to another nucleus of any given type (A), and the proportions of nuclei of different types that immediately surround a Type-2 nucleus (B). Here the analysis was performed on the entire dataset pooled from 14 experimental fields; see Figure S2 for results performed on individual embryos as in panel (B).
To further confirm the clustering property of Type-2 nuclei observed by Ripley’s K-function, we performed another analysis. If Type-2 nuclei in our experimental fields do form clusters that are not random in nature, we would expect to observe a difference between our experimental data (blue in Figures 5B and S2) and randomly simulated data with regard to the distributions of the types of nuclei immediately surrounding a Type-2 nucleus. We calculated the envelopes of the 95% confidence interval (CI) from either a theoretical model (red in Figures 5B and S2) or random simulations (green in Figures 5B and S2). The theoretical model and random simulations are detailed in Materials and Methods. Among nuclei surrounding a Type-2 nucleus, the proportions of both Type-0 and Type-2 nuclei lie outside the 95% CI, with the former being smaller and latter being greater than random. Equivalent results were also obtained for the proportions of different types of nuclei that surround a Type-0 nucleus, but not when evaluating the neighborhood of a Type-1 nucleus (data not shown). Together these findings suggest that, in early cycle 14 embryos, transcriptional states of hb genes are likely to be correlated among immediate-neighboring nuclei.
A Nuclear Lineage Model Evaluating the Outcome of Inheritance of Transcriptional States
A simplified, idealized model of inheritance or memorization can be stated as follows, although the actual biological system is likely much more complex. A nucleus derived from a parental nucleus upon a mitotic cycle adopts both the parental location (i.e., either at or immediately adjacent to this location) and the parental state of transcription (i.e., of the same type). To determine whether this simple model may be able to capture some of the key features of the experimentally observed spatial distributions of the nuclei, we performed simulations using a series of nuclear lineage patterns from cycle 11 to cycle 14. Using the positions of experimentally observed cycle-14 nuclei (the mold) shown in Figure 1B, we assigned randomly all nuclei into 8-neighbor lineage clusters (see Materials and Methods and also Figure S3). Based on experimental observations , we used a 90% transcribing probability for hb gene copies at cycles 11 and 12, and 50% at cycles 13 and 14. This decrease in transcribing probability was shown to coincide with a decay of the maternal Hb, marking a transition from a virtually uniform (also referred to as “synchronous”) to stochastic hb transcription . The number distributions of the experimentally observed three types of nuclei within the hb transcribing domain are well described by a binomial function irrespective of the overall transcribing probability , .
We performed 99 simulations either with or without inheritance. When inheritance of transcriptional states between mitotic cycles was not implemented, all hb gene copies were allowed to undergo active transcription in a random and independent manner. In this case, the number of active copies in each nucleus followed a binomial distribution, which was a function of transcribing probability at each cycle (Figure 6A). When inheritance was implemented, the two daughter nuclei (at each new cycle) were assigned to exist as the same type as their parental nucleus. Active copies were randomly switched off to achieve the experimentally observed reduction in transcribing probability at cycle 13 (Figure 6B). To evaluate the spatial distributions of the simulated data, we computed the H-function profiles. Without inheritance, Types-0 and -2 nuclei behave closely to CSR (Figure 6C; see also Figure S4A for one simulated pattern). With inheritance, the H-functions of these two types of nuclei become positive values within the r range of 6~78 µm. These simulated properties are robust to the choices of transcribing probabilities: reducing the overall transcribing probability at cycle 11 and 12 from 90% to 60% or varying the probability at cycles 13 and 14 between 40% and 60% did not alter the outcomes (data not shown). We note that the simulated clustering (Figure 6D; see also Figure S4B) is tighter than observed in the experimental data (Figure 3A and 3C), with both greater H(r) values and a narrower r range, especially for Type-2 nuclei. This difference suggests that the actual biological system may not follow the strict inheritance as implemented in our idealized model (see below for further discussions). This difference may also be contributed by additional underlying differences between the actual biological system and our simulated system. For example, the 8-nuclear lineage clusters used in our simulations (3.4±0.5% EL; see also Figure S3A) are tighter than the experimentally observed clusters (4.8±1% EL, ).
Shown are simulated examples (A and B) and the H-function profiles (C and D) of the nuclear types with (B and D) or without (A and C) inheritance of transcriptional states. These depicted examples are for illustrative purposes by showing the simulated transcription state lineages with nuclei aligned in one dimension; see Figure S3 for simulated examples of 8-nuclei lineage clusters arranged in two-dimensional fields. For the initial conditions at cycle 11, the probability of hb busting is set at a high level as observed experimentally (see also Materials and Methods). With an overall transcribing probability of 90%, the majority of nuclei are Type 2 (81%), with Type-1 being 18% and Type-0 being 1% (not shown in the diagrams that are for illustrative purposes only). At cycle 13, a binomial process is implemented in our simulation that randomly shuts off active hb copies so that the overall transcribing probability becomes 50% (with Types-0, -1 and -2 being 25%, 50% and 25%, respectively). The type of each nucleus in (A) is assigned by the binomial random number generator of MATLAB.
Understanding how the spatial information is encoded by morphogen gradients and decoded by cells is of fundamental importance to developmental biology and biophysics –. Recent studies, fueled by both technological and conceptual advances, have uncovered important insights into this problem. Most of the existing studies have focused primarily and directly on how the boundaries of gene expression patterns are established. Unlike these studies, our current study focuses on a distinct, but related (see below), aspect of transcriptional responses to morphogen action. Specifically, we investigate the spatial distributions of nuclei exclusively within the hb transcribing domain of Drosophila embryos. Our previous studies of hb transcriptional states were based on grouping of nuclei into bins along the AP axis to permit direct evaluations of the threshold response to the Bcd gradient on an embryonic scale , . A bin of 2% EL contains on average <2 nuclei along the AP axis. Meanwhile, the nuclear lineage from cycle 11 to 14 has an average width of 4.8±1% EL . Thus the binning would arbitrarily break the territories of any nuclear lineages and, consequently, lead to a loss of the full spatial information on a finer scale. A bin within in the transcribing domain can cross >10 different territories , effectively eliminating (through averaging within the bin) any nucleus-to-nucleus correlations that may exist within individual territories. This feature is suitable for specifically evaluating transcriptional states along the AP axis. In contrast to the binning method, our current study preserves the native spatial information about the nuclei in the monolayer within the hb transcribing domain. This permits the extraction of the features about the spatial relationships among the nuclei with regard to their hb transcriptional states.
Our experimental results show that, due to the nearly uniform size of the nucleus and the tight packing against one another in the blastoderm embryo, the nuclei themselves are regularly spaced when evaluated on a scale at or below the mean nuclear diameter (Figure 2). However, on larger scales, nuclei with both or neither hb copies undergoing active transcription (i.e., Type-2 and Type-0) exhibit properties of spatial clustering. These results can be explained by a model where transcriptional states of the nuclei are inherited between mitotic cycles (Figure 6D). Thus, although our current analysis focuses exclusively on properties within the transcribing domain, it illuminates on the potential inheritance of transcriptional states, an issue fundamental to how gene expression boundaries are formed in response to the Bcd concentration gradient during development (see Introduction). Our findings are supportive of a hypothesis that active hb transcription may not depend on de novo sensing–at all mitotic cycles–of the nuclear Bcd concentration to form an expression boundary .
In our analysis of the experimental data, nuclei with a single active hb copy showed little or no evidence of spatial clustering. We currently do not fully understand why Type-1 nuclei behave differently from Type-0 and Type-2 nuclei. We note that our simulations based on the experimentally obtained molds also revealed a lack of spatial clustering of Type-1 nuclei (Figure 6D). Type-1 and Type-2 nuclei are similar in the sense that they both have at least one transcribing hb copy and are thus transcriptionally “active” , . But they are vastly different from a statistical point of view. The binomial distribution of the observed probability of the three types of nuclei within the transcribing domain dictates that Type-1 nuclei are the largest fraction at cycle 14 (when the overall transcribing probability is 50%; see Figure 6 legend for expected fractions of the three types based on binomial distribution). Thus, among all three types, Type-1 nuclei have the highest probability of being involved in potential type switches. These switches would likely take place at random locations in the embryo, which may further contribute to both the observed lack of clustering of Type-1 nuclei and the relatively loose clustering of the other two types.
We note that, since our statistical analyses cannot draw direct conclusions about causal relationships, experimental validations of inheritance of the transcriptional states for Bcd target genes such as hb will await technological advances in the future. We also note that there are documented precedents for inheritance of transcriptional states between mitotic cycles. For example, studies of the transcriptional states of the globin genes within clonally-expanded populations of mouse cells suggested that stochastic decisions can be inherited and maintained between mitotic cycles . Positional effect variegation (PEV) is another example of mitotic inheritance where the chromatin structure is known to play a role , . In neither case, to our knowledge, has the timing of the onset of active transcription of individual gene copies upon exiting the mitosis been carefully evaluated. While these and other examples of inheritance may be related to the inheritance that we investigate here, the timing of the onset of active transcription after mitosis in our case is very quick, which is on scales that are likely measured by seconds or tens of seconds as opposed to minutes or tens of minutes , , . It has been proposed that promoters with stalled RNAP can respond to induction more quickly , . There is evidence that the hb P2 promoter, similar to promoters of many other patterning genes , contains stalled RNAP . If the stalled RNAP were responsible for the inheritance that we investigate here, how would it arrive at the gene promoter so quickly (after mitosis) in the first place? Stalled RNAP probably could not have survived the mitotic process. It has been proposed that transcription takes place at discrete nuclear locations referred to as the transcription factories , . Could the actively-transcribing hb copies, i.e., the intron dots, take residence in these factories? Could the Bcd molecules that are enriched at the intron dot locations , and/or maternal Hb which plays a role in increasing the transcribing probability of hb in earlier cycles , also take residence in these factories? It remains a significant future challenge to understand whether and in what fashion these factories might contribute, if at all, to the memorization or inheritance of transcriptional states between mitotic cycles. Despite unresolved questions, our current study represents an important step toward understanding both the decoding of positional information during development and the fundamental transcription process.
Location of the experimental field within the contours of all 14 embryos. The dots in different colors represent the nuclear positions at the outer edges of 14 embryos; the solid box represents the field, where our analyses were performed.
Estimation of proportions of nuclei surrounding a Type-2 nucleus for individual embryos. Same as Figure 5B but for individual embryos.
Simulated nuclear lineages. Shown are nuclear lineages generated from one random simulation based on the mold shown in Figure 1B. We performed 99 simulations in total. In each panel, neighboring nuclear lineages are distinguished by different colors. (A-C) represent 8-nuclei lineages from cycle 11 to cycle 14, 4-nuclei lineages from cycle 12 to cycle 14, and 2-nuclei lineages from cycle 13 to cycle 14, respectively. See Materials and Methods for details.
Simulated nuclear patterns of transcriptional states. Shown are examples of simulated fields exhibiting the transcriptional states of the nuclei, with (B) or without (A) inheritance. These two simulations were based on the exact nuclear lineage assignments shown in Figure S1.
Nuclear numbers of each type in the experimental fields of all 14 embryos. For each type of nuclei, the mean and the standard deviation among 14 embryos are given in the last column.
We thank Jie Ren for generating the experimental results and confocal images used in the current analysis, Junbo Liu for sharing unpublished information, and Jingyuan Deng and other members of our CCHMC groups for discussions and technical assistance. We thank the anonymous reviewer and the editor for their constructive suggestions to improve the presentation of this work.
Conceived and designed the experiments: FH JM. Performed the experiments: FH. Analyzed the data: FH. Contributed reagents/materials/analysis tools: FH. Wrote the paper: FH JM.
- 1. Weake VM, Workman JL (2010) Inducible gene expression: diverse regulatory mechanisms. Nat Rev Genet 11: 426–437. doi: 10.1038/nrg2781
- 2. Ma J (2011) Transcriptional activators and activation mechanisms. Protein Cell 2: 879–888. doi: 10.1007/s13238-011-1101-7
- 3. Nechaev S, Adelman K (2011) Pol II waiting in the starting gates: Regulating the transition from transcription initiation into productive elongation. Biochim Biophys Acta 1809: 34–45. doi: 10.1016/j.bbagrm.2010.11.001
- 4. Boettiger AN, Levine M (2009) Synchronous and stochastic patterns of gene activation in the Drosophila embryo. Science 325: 471–473. doi: 10.1126/science.1173976
- 5. He F, Ren J, Wang W, Ma J (2011) A multiscale investigation of bicoid-dependent transcriptional events in Drosophila embryos. PLoS One 6: e19122. doi: 10.1371/journal.pone.0019122
- 6. To TL, Maheshri N (2010) Noise can induce bimodality in positive transcriptional feedback loops without bistability. Science 327: 1142–1145. doi: 10.1126/science.1178962
- 7. Zobeck KL, Buckley MS, Zipfel WR, Lis JT (2010) Recruitment timing and dynamics of transcription factors at the Hsp70 loci in living cells. Mol Cell 40: 965–975. doi: 10.1016/j.molcel.2010.11.022
- 8. Stavreva DA, Varticovski L, Hager GL (2012) Complex dynamics of transcription regulation. Biochim Biophys Acta 1819: 657–666. doi: 10.1016/j.bbagrm.2012.03.004
- 9. He F, Ren J, Wang W, Ma J (2012) Evaluating the Drosophila Bicoid morphogen gradient system through dissecting the noise in transcriptional bursts. Bioinformatics 28: 970–975. doi: 10.1093/bioinformatics/bts068
- 10. Eldar A, Elowitz MB (2010) Functional roles for noise in genetic circuits. Nature 467: 167–173. doi: 10.1038/nature09326
- 11. Tkacik G, Gregor T, Bialek W (2008) The role of input noise in transcriptional regulation. PLoS One 3: e2774. doi: 10.1371/journal.pone.0002774
- 12. Chubb JR, Liverpool TB (2010) Bursts and pulses: insights from single cell studies into transcriptional mechanisms. Curr Opin Genet Dev 20: 478–484. doi: 10.1016/j.gde.2010.06.009
- 13. Shermoen AW, O'Farrell PH (1991) Progression of the cell cycle through mitosis leads to abortion of nascent transcripts. Cell 67: 303–310. doi: 10.1016/0092-8674(91)90182-x
- 14. Sutherland H, Bickmore WA (2009) Transcription factories: gene expression in unions? Nat Rev Genet 10: 457–466. doi: 10.1038/nrg2592
- 15. Bothma JP, Magliocco J, Levine M (2011) The snail repressor inhibits release, not elongation, of paused Pol II in the Drosophila embryo. Curr Biol 21: 1571–1577. doi: 10.1016/j.cub.2011.08.019
- 16. McHale P, Mizutani CM, Kosman D, MacKay DL, Belu M, et al. (2011) Gene length may contribute to graded transcriptional responses in the Drosophila embryo. Dev Biol 360: 230–240. doi: 10.1016/j.ydbio.2011.08.016
- 17. Porcher A, Abu-Arish A, Huart S, Roelens B, Fradin C, et al. (2010) The time to measure positional information: maternal hunchback is required for the synchrony of the Bicoid transcriptional response at the onset of zygotic transcription. Development 137: 2795–2804. doi: 10.1242/dev.051300
- 18. Liu J, Ma J (2013) Uncovering a dynamic feature of the transcriptional regulatory network for anterior-posterior patterning in the Drosophila embryo. PLoS One (in revision).
- 19. Driever W, Nusslein-Volhard C (1988) The bicoid protein determines position in the Drosophila embryo in a concentration-dependent manner. Cell 54: 95–104. doi: 10.1016/0092-8674(88)90183-3
- 20. Struhl G, Struhl K, Macdonald PM (1989) The gradient morphogen bicoid is a concentration-dependent transcriptional activator. Cell 57: 1259–1273. doi: 10.1016/0092-8674(89)90062-7
- 21. Ephrussi A, St Johnston D (2004) Seeing is believing: the bicoid morphogen gradient matures. Cell 116: 143–152. doi: 10.1016/s0092-8674(04)00037-6
- 22. Jaeger J, Surkova S, Blagov M, Janssens H, Kosman D, et al. (2004) Dynamic control of positional information in the early Drosophila embryo. Nature 430: 368–371. doi: 10.1038/nature02678
- 23. Holloway DM, Lopes FJ, da Fontoura Costa L, Travencolo BA, Golyandina N, et al. (2011) Gene expression noise in spatial patterning: hunchback promoter structure affects noise amplitude and distribution in Drosophila segmentation. PLoS Comput Biol 7: e1001069. doi: 10.1371/journal.pcbi.1001069
- 24. de Lachapelle AM, Bergmann S (2010) Precision and scaling in morphogen gradient read-out. Mol Syst Biol 6: 351. doi: 10.1038/msb.2010.7
- 25. White MA, Parker DS, Barolo S, Cohen BA (2012) A model of spatially restricted transcription in opposing gradients of activators and repressors. Mol Syst Biol 8: 614. doi: 10.1038/msb.2012.48
- 26. Zhang L, Radtke K, Zheng L, Cai AQ, Schilling TF, et al. (2012) Noise drives sharpening of gene expression boundaries in the zebrafish hindbrain. Mol Syst Biol 8: 613. doi: 10.1038/msb.2012.45
- 27. Tostevin F, ten Wolde PR, Howard M (2007) Fundamental limits to position determination by concentration gradients. PLoS Comput Biol 3: e78. doi: 10.1371/journal.pcbi.0030078.eor
- 28. Bialek W, Setayeshgar S (2005) Physical limits to biochemical signaling. Proc Natl Acad Sci U S A 102: 10040–10045. doi: 10.1073/pnas.0504321102
- 29. Gregor T, Tank DW, Wieschaus EF, Bialek W (2007) Probing the limits to positional information. Cell 130: 153–164. doi: 10.1016/j.cell.2007.05.025
- 30. He F, Saunders TE, Wen Y, Cheung D, Jiao R, et al. (2010) Shaping a morphogen gradient for positional precision. Biophys J 99: 697–707. doi: 10.1016/j.bpj.2010.04.073
- 31. Ripley BD (1977) Modeling Spatial Patterns. Journal of the Royal Statistical Society Series B-Methodological 39: 172–212.
- 32. Diggle PJ (1986) Displaced Amacrine Cells in the Retina of a Rabbit - Analysis of Bivariate Spatial Point Pattern. Journal of Neuroscience Methods 18: 115–125. doi: 10.1016/0165-0270(86)90115-9
- 33. Luengo Hendriks CL, Keranen SV, Fowlkes CC, Simirenko L, Weber GH, et al. (2006) Three-dimensional morphology and gene expression in the Drosophila blastoderm at cellular resolution I: data acquisition pipeline. Genome Biol 7: R123.
- 34. Clark PJ, Evans FC (1954) Distance to nearest neighbor as a measure of spatial relationships in populations. Ecology 35: 445–453. doi: 10.2307/1931034
- 35. Okabe A, Boots B, Sugihara K, Chiu SN (2000) Spatial tessellations: concepts and applications of Voronoi diagrams. New York: Wiley.
- 36. Burz DS, Rivera-Pomar R, Jackle H, Hanes SD (1998) Cooperative DNA-binding by Bicoid provides a mechanism for threshold-dependent gene activation in the Drosophila embryo. EMBO J 17: 5998–6009. doi: 10.1093/emboj/17.20.5998
- 37. Ma X, Yuan D, Diepold K, Scarborough T, Ma J (1996) The Drosophila morphogenetic protein Bicoid binds DNA cooperatively. Development 122: 1195–1206.
- 38. Perry MW, Boettiger AN, Levine M (2011) Multiple enhancers ensure precision of gap gene-expression patterns in the Drosophila embryo. Proc Natl Acad Sci U S A 108: 13570–13575. doi: 10.1073/pnas.1109873108
- 39. Blankenship JT, Wieschaus E (2001) Two new roles for the Drosophila AP patterning system in early morphogenesis. Development 128: 5129–5138.
- 40. Miles CM, Lott SE, Hendriks CL, Ludwig MZ, Manu, et al (2011) Artificial selection on egg size perturbs early pattern formation in Drosophila melanogaster. Evolution 65: 33–42. doi: 10.1111/j.1558-5646.2010.01088.x
- 41. Raven MA, Stagg SB, Reese BE (2005) Regularity and packing of the horizontal cell mosaic in different strains of mice. Vis Neurosci 22: 461–468. doi: 10.1017/s0952523805224070
- 42. Eglen SJ, Raven MA, Tamrazian E, Reese BE (2003) Dopaminergic amacrine cells in the inner nuclear layer and ganglion cell layer comprise a single functional retinal mosaic. J Comp Neurol 466: 343–355. doi: 10.1002/cne.10891
- 43. Kiskowski MA, Hancock JF, Kenworthy AK (2009) On the use of Ripley's K-function and its derivatives to analyze domain size. Biophys J 97: 1095–1103. doi: 10.1016/j.bpj.2009.05.039
- 44. Bell ML, Earl JB, Britt SG (2007) Two types of Drosophila R7 photoreceptor cells are arranged randomly: a model for stochastic cell-fate determination. J Comp Neurol 502: 75–85. doi: 10.1002/cne.21298
- 45. Arias AM, Hayward P (2006) Filtering transcriptional noise during development: concepts and mechanisms. Nat Rev Genet 7: 34–44. doi: 10.1038/nrg1750
- 46. Kerszberg M, Wolpert L (2007) Specifying positional information in the embryo: looking beyond morphogens. Cell 130: 205–209. doi: 10.1016/j.cell.2007.06.038
- 47. Lewis J (2008) From signals to patterns: space, time, and mathematics in developmental biology. Science 322: 399–403. doi: 10.1126/science.1166154
- 48. Lander AD (2011) Pattern, growth, and control. Cell 144: 955–969. doi: 10.1016/j.cell.2011.03.009
- 49. de Krom M, van de Corput M, von Lindern M, Grosveld F, Strouboulis J (2002) Stochastic patterns in globin gene expression are established prior to transcriptional activation and are clonally inherited. Mol Cell 9: 1319–1326. doi: 10.1016/s1097-2765(02)00558-0
- 50. Milot E, Strouboulis J, Trimborn T, Wijgerde M, de Boer E, et al. (1996) Heterochromatin effects on the frequency and duration of LCR-mediated gene transcription. Cell 87: 105–114. doi: 10.1016/s0092-8674(00)81327-6
- 51. Le HD, Donaldson KM, Cook KR, Karpen GH (2004) A high proportion of genes involved in position effect variegation also affect chromosome inheritance. Chromosoma 112: 269–276. doi: 10.1007/s00412-003-0272-2
- 52. Levine M (2011) Paused RNA polymerase II as a developmental checkpoint. Cell 145: 502–511. doi: 10.1016/j.cell.2011.04.021
- 53. Zeitlinger J, Stark A, Kellis M, Hong JW, Nechaev S, et al. (2007) RNA polymerase stalling at developmental control genes in the Drosophila melanogaster embryo. Nat Genet 39: 1512–1516. doi: 10.1038/ng.2007.26
- 54. Li H, Rodriguez J, Yoo Y, Shareef MM, Badugu R, et al. (2011) Cooperative and antagonistic contributions of two heterochromatin proteins to transcriptional regulation of the Drosophila sex determination decision. PLoS Genet 7: e1002122. doi: 10.1371/journal.pgen.1002122
- 55. Osborne CS, Chakalova L, Brown KE, Carter D, Horton A, et al. (2004) Active genes dynamically colocalize to shared sites of ongoing transcription. Nat Genet 36: 1065–1071. doi: 10.1038/ng1423
- 56. Razin SV, Gavrilov AA, Pichugin A, Lipinski M, Iarovaia OV, et al. (2011) Transcription factories in the context of the nuclear and genome organization. Nucleic Acids Res 39: 9085–9092. doi: 10.1093/nar/gkr683