Dickkopf1 - A New Player in Modelling the Wnt Pathway

The Wnt signaling pathway transducing the stabilization of β-catenin is essential for metazoan embryo development and is misregulated in many diseases such as cancers. In recent years models have been proposed for the Wnt signaling pathway during the segmentation process in developing embryos. Many of these include negative feedback loops where Axin2 plays a key role. However, Axin2 null mice show no segmentation phenotype. We therefore propose a new model where the negative feedback involves Dkk1 rather than Axin2. We show that this model can exhibit the same type of oscillations as the previous models with Axin2 and as observed in experiments. We show that a spatial Wnt gradient can consistently convert this temporal periodicity into the spatial periodicity of somites, provided the oscillations in new cells arising in the presomitic mesoderm are synchronized with the oscillations of older cells. We further investigate the hypothesis that a change in the Wnt level in the tail bud during the later stages of somitogenesis can lengthen the time period of the oscillations and hence the size and separation of the later somites.


Introduction
A segmented body plan is a fundamental characteristic feature of vertebrates.The process of segmentation is carried out by a combination of changes in gene expression and relative anterior-posterior cell position in the presomitic mesoderm (PSM) (?).In the anterior end of the embryo the somites are segmented at equally separated time points with species dependent periods.In mice the period is around 120 min and in frogs it is around 90 min.
In (?) Cooke and Zeeman proposed the clock and wavefront model to describe the segmentation process.The idea is that locally coupled oscillators are controlled by a morphogen gradient in the PSM.The oscillators are the clocks providing temporal information, e.g., cycle state, and the morphogen gradient is the wavefront providing spatial information about axial position.Until now three major pathways controlling the segmentation process have been found: the Notch, Wnt and FGF pathways.They all have target genes, which oscillates and, interestingly, Wnt target genes oscillate out of phase with Notch and FGF target genes (?).These three pathways could be the clocks.There is a decreasing gradient of wnt3a and fibroblast growth factor 8 (fgf8) starting from the tail bud through the PSM (??).Fgf8 seems to be regulated by Wnt (?), thus Wnt could form the main wavefront.The actual setting of the somites happens at the determination front, where the gradient drops below a certain threshold.Cells past the determination front become permissive to form somites depending on their phase of oscillation (??).
In 2003 it was discovered by Aulehla et al. that Axin2 oscillates during the segmentation process in developing mouse embryos.Since their discovery several models for the Wnt oscillator have been proposed (???) with Axin2 as the main variable.We propose a new model for the Wnt pathway with Dickkopf1 (Dkk1) as the core of a negative feedback loop.The reason for changing from Axin2 to Dkk1 is that Axin2 deleted mice (Axin2 − / − ) have no dramatic phenotype (?).Dkk1 has an oscillatory behavior during the segmentation process in mouse embryos (?) and lack of Dkk1 in mouse embryos results in anterior head truncations and neonatal death (?).In mice, for which we want to model the segmentation process, several Dickkorpf genes are known.Dkk1 is the most widely studied of these.

Results and Discussion
Modeling the Wnt/β-catenin pathway During Wnt signaling β-catenin interacts with the TCF/LEF-1 DNA-binding proteins to promote transcription of Wnt target genes (?).As for Axin2 the transcription factor for Dkk1 is β-catenin (??).After transcription and translation Dkk1 goes through the cellular membrane where it can bind to the extracellular domains of the low-density lipoprotein receptor-related protein 5 and 6 (LRP5/6).When bound to LRP5/6, Dkk1 acts as an inhibitor of Wnt signaling by blocking the association between Wnt, Frizzled (Fz) and LRP5/6 (?).Wnt acts as an inducer for the formation of this complex and Dkk1 is a competitor to this induction (??).
It has been proposed that the Wnt signal is transduced through the cell membrane by the binding of Dishevelled (Dsh) to the intracellular domain of the Fz receptor(?).Axin and Dsh can bind together via their DIX domains (?) and they co-localize at the membrane (?) during Wnt signaling.Therefore Dsh bound to Fz may recruit Axin bound to GSK3β to the LRP5/6 receptor (?), where a phosphorylation of LRP5/6 is initiated.The LRP5/6 receptor has a binding site for Axin and upon Wnt signalling GSK3β (bound to Axin) phosphorylates LRP5/6, which requires Axin (?).The phosphorylated LRP5/6 receptor may be able to recruit and more efficiently bind the Axin-GSK3β complex to the membrane and the phosphorylation process is thereby amplified (?).
At the cell membrane Axin is phosphorylated by GSK3β and then degraded (??).The degradation of Axin leads to a decrease in the formation of the destruction complex comprised of β-catenin, the kinases glycogen synthase kinase 3 (GSK3β), casein kinase (CKIα), the scaffolding proteins Axin and adenomatous polyposis coli (APC).In the destruction complex β-catenin gets phosphorylated and subsequently degraded.
Interestingly enough GSK3β plays a dual role in controlling the Wnt signal.When the Wnt signal is off then GSK3β phosphorylates β-catenin in the destruction complex and when the signal is on then it phosphorylates Axin at the LRP5/6 receptor.Whether it is the same or distinct Axin-GSK3β complexes that carry out the phosphorylation of β-catenin and Axin is unknown (?).
The dynamics of CK1α and APC are included in the parameters governing the destruction complex, the TCF/LEF-1 dynamics are contained in the transcription of Dkk1, and the dynamics of Dsh are included in the formation of the complex consisting of Axin, GSK3β and LRP5/6 at the cell membrane.Figure 1 shows a simplified diagram of the proposed model and the associated equations.C, [GA], G, B, L, D m , D, [LD], A and [LGA] are the concentrations of the destruction complex, GSK3β-Axin complex, GSK3β, β-catenin, LRP5/6, Dkk1 mRNA, Dkk1 protein, Dkk1-LRP5/6 complex, Axin and LRP5/6-Axin-GSK3β complex.The formation and breaking of a complex X are denoted by c f X and c bX , respectively.The transcription and translation rate of Dkk1 are given by the parameters c tsc and c tsl , respectively.The Hill coefficient on B regarding the transcription of D m is associated to the amount of cooperativity between β-catenin and the TCF/LEF-1 complex.For example, no cooperativity would result in a Hill coefficient of one.
The concentration of GSK3β has been shown to be extremely stable (?) and consequently its total concentration, GSK3β tot , is assumed to be constant during the involved time scales.The same assumption goes for the total concentration of LRP5/6, since the half-life of LRP6 is around 4.7 hours (?).Therefore we have not included any source or sink for the concentrations of G and L. Only constitutive sources of β-catenin and Axin, respectively, are included in the model, because free (unphosphorylated) β-catenin and Axin is stable (?).

Determination of the parameter values
In ? a model for the Wnt pathway in Xenopus is proposed and from this article estimates for the parameters , α, GSK3β tot , S A and S B are taken.Dkk1 binds to LRP5/6 with a high affinity, since the dissociation constant has been measured to be around K [LD] = 0.4 − 0.5 nM (??).
The other parameters are estimated to produce oscillations with a period of around 120 min and a very low concentration of total Axin as found by ?.The low concentration of total Axin is thought to act as a buffer to changes in the concentration of the other constituents, which may also take part in other signalling pathways (?).
In the activation of transcription of Dkk1 by β-catenin we assume cooperativity between β-catenin and TCF/LEF.The Hill coefficient is set to three.A Hill coefficient of one and two returns no sustained oscillations, when the model is simulated.In other models of the Wnt pathway Hill coefficients of values 2 and 5 have been applied (??).

Oscillations of the mRNA, protein and complex levels
The parameter values produce oscillations of the involved constituents with a period of around 120 min as seen in Fig. 1.In addition, the concentration of Axin is very low compared to the concentrations of the other variables, as discussed above.The different steps in the proposed Wnt model comply with the phase shifts in the oscillations of the different components.An increase in the [LGA] concentration leads to a decrease in the Axin concentration.This decrease causes a reduction in the concentration of the destruction complex and consequently an increase in the β-catenin concentration.This increase will after a while cause an increase in the Dkk1 concentration.The concentrations of GSK3β and [GA] are mirrors of each other, since a high concentration of will cause [GA] to leave less free GSK3β behind.
? found no significant oscillations in the level of β-catenin.For our choice of parameters β-catenin show oscillatory behavior with an amplitude of approximately 5 nM, which is not significantly low.A different set of parameters could possibly give a smaller amplitude of β-catenin but the general results, presented later, is not significantly altered by this.

A spatial Wnt gradient induced by time variation
The described model could be the clock of the clock and wavefront model as mentioned above.The wavefront we introduce as a time dependent Wnt gradient through a time dependence of the parameter c f [LGA] .The reason for making c f [LGA] time dependent is clear if we introduce the variable [LW ] describing the binding of Wnt (W ) to the LRP5/6 receptor together with Frizzled.The temporal change of [LW ] is described as and if we assume steady state for the binding of Wnt to LRP5/6 then Inserting this L estimate into the concentration of [LGA], it can be seen that the formation of [LGA] will be time dependent, with To include the decreasing gradient of Wnt observed in mouse embryos we use a decreasing value of c f [LGA] .If we assume that Wnt primarily decreases through a diffusive process, then the profile of c f [LGA] will be Gaussian.In Fig. 3 the amplitude and period of Dkk1 concentration is plotted as c f [LGA] changes.The green line refers to the reference value of c f [LGA] ; reducing it results in smaller amplitudes.The period seems to be almost unaltered, which complies with the findings of constant periods of the oscillations in the PSM.
The length of the PSM is approximately constant during the formation of the first 15-20 somites in mice embryos (?).The same is almost true for the size of the somites.From ? one can find that a cell budded off in the tail bud will be segmented within the order of ∼ 1100 − 1300 mins.
It has been measured that FGF (regulated by Wnt) exhibits a gradient in the PSM with a fold change of 2-5 (?).The fold change used for c f [LGA] is 2. Assuming that Wnt diffuses through the PSM and setting the final value of c f [LGA] equal to 2.5 1 nMmin enable us to calculate a Gaussian profile of c f [LGA] representing the Wnt gradient.In Supplementary Fig. 1, Gaussian profiles of c f [LGA] are plotted with different initial values.A decreasing value of c f [LGA] in the PSM will give rise to smaller amplitudes and slightly longer periods.

Synchronization of neighbouring cells
The elongation of the embryo is modelled by simulating cells budding off from the tail bud at regular time intervals, R. Elongation is considered in only one dimension.In Fig. 2 time series for cells with R = 10 min are shown.The final level of Dkk1 at the determination front is oscillating (see Fig. 2A) with a period of 120 min.However, these oscillations are only possible if the onset of oscillations in the cells are synchronized, i.e. if the state for cell i at time t is denoted S i (t), then S i+1 (0) = S i (R).This means that the initial state of a cell budded off from the tail bud equals that of its anterior neighbour.If all the cells have the same initial state then the levels of Dkk1 at the determination front will also be equal, see Fig. 2B.
One could imagine that the synchronization of the clocks in each cell is not perfect.If the initial state of cell i + 1 is chosen randomly as )), then the period is almost unaltered (see Supplementary Fig. 2A).Here U(−x, x) denotes random numbers uniformly distributed between −x and x.If the initial state is randomly chosen within the range of Dkk1 levels, then the oscillations of Dkk1 at the determination front are disrupted and no periodicity is visual (see Supplementary Fig. 2B).Thus, the system seems to be robust to small changes in the synchronization.
Diffusion of Wnt is thought to be the main reason for the Wnt gradient in the PSM.Thus a Gaussian profile of c f [LGA] is used.The profile of c f [LGA] would still be Gaussian if we also included a half-life of Wnt.If the Wnt gradient was only controlled by the half-life of Wnt then the levels of Dkk1 still oscillate at the determination front with a period of 120 min, when the model is simulated as above with R = 10 min.

The period of somite formation increases
It is known that the period of somite formation increases as the embryo elongates (?).This led ? to propose that the level of Wnt production in the tail bud increases with time.To test their hypothesis we let the initial value of c f [LGA] increase with time keeping the level at the determination front fixed at 2.5 1 nMmin as above.Simulating an increase in the Wnt level of the tail bud results in longer periods of the Dkk1 levels at the determination front, see Supplementary Fig. 3. Thus, we can confirm the hypothesis of ?.Note that a steeper increase of the Wnt level in the tail bud results in a greater extension of the periods (compare Supplementary Figs.3A and 3B).

Outlook
We propose a Wnt model with Dkk1 as the core for the negative feedback loop, which exhibits sustained oscillations of Dkk1.By simulation of the elongating embryo we were able to test the importance of synchronization between neighbouring cells.In addition, we could also show that a small change in the synchronization did not significantly disrupt smooth oscillations of Dkk1 levels at the determination front.The clock and wavefront model was investigated through a decreasing Wnt gradient in the PSM and we could confirm the hypothesis that the Wnt level in the tail bud increases over time resulting in longer periods of the oscillations.Further aspects of the growing embryo could be incorporated into the model, like a varying PSM size and a rate of cells budding off from the tail bud.Included are only members of the Wnt pathway that are important for the understanding of the negative feedback loop.When the Wnt signal is on, then Axin gets degraded at the LRP5/6 complex and β-catenin can act as a transcription factor of the Wnt inhibitor Dkk1.Vice-versa, when the Wnt signal is off, due to inhibition by Dkk1, then β-catenin gets degraded.Top,right Simulated time series for a selection of variables from the model listed in the bottom panel.The total level of Axin (magenta) is low, which complies with the findings of ?. Bottom: The terms describing the time dependence of the variables.period= 130 min

Figure 1 :
Figure 1: Top,left: A diagram of the Wnt model with a feed-back loop over Dkk1.Included are only members of the Wnt pathway that are important for the understanding of the negative feedback loop.When the Wnt signal is on, then Axin gets degraded at the LRP5/6 complex and β-catenin can act as a transcription factor of the Wnt inhibitor Dkk1.Vice-versa, when the Wnt signal is off, due to inhibition by Dkk1, then β-catenin gets degraded.Top,right Simulated time series for a selection of variables from the model listed in the bottom panel.The total level of Axin (magenta) is low, which complies with the findings of ?. Bottom: The terms describing the time dependence of the variables.

Figure 2 :
Figure 2: Top: Time series for the Dkk1 concentration.Space is introduced by letting a cell bud off from the tail bud every 10th min.Thus the cells move relatively in the PSM.At the determination front the oscillations arrest.Bottom: The level of Dkk1 at the determination front.(A) R = 10 min.(B) R = 10 min.(C) R = 10 min and the cells have the same initial level of Dkk1 TB: tail bud.S i : Somite i, where S 0 is the newly formed somite.

Figure 3 :
Figure 3: The parameter c f [LGA] is changed, and the amplitude (solid) and the period (dashed) are calculated.The green vertical line denotes the value of c f [LGA] in the reference state.