OCT4 expression in human embryonic stem cells: spatio-temporal dynamics and fate transitions
Laura E Wadkin, Sirio Orozco-Fuentes, Irina Neganova, Majlinda Lako, R A Barrio, A W Baggaley, Nicholas G Parker, Anvar Shukurov
OOCT4 expression in humanembryonic stem cells: temporalproperties and regulatorydynamics
L E Wadkin , S Orozco-Fuentes , , INeganova , M Lako , N G Parker and AShukurov School of Mathematics, Statistics and Physics, NewcastleUniversity, Newcastle upon Tyne, UK. Department of Mathematics, Physics and ElectricalEngineering, Northumbria University, Newcastle upon Tyne,UK. Institute of Cytology, RAS St Petersburg, Russia Biosciences Institute, Faculty of Medical Sciences, NewcastleUniversity, Newcastle upon Tyne, UK.E-mail: [email protected]
May 2020
Abstract.
Pluripotency is the defining characteristic of humanembryonic stem cells (hESCs) allowing them to differentiateinto any somatic cell in the human body. For thepromising clinical applications of hESCs, improved regulationof pluripotency and differentiation trajectories of their coloniesis required. The pluripotency transcription factors (PTFs)which regulate pluripotency are inherently stochastic (with smallfluctuations impacting cell fate), inherited asymmetrically uponcell divisions, and more similar in closely related (separated byless division events) cells. Here we use time-lapse experimentaldata of OCT4 fluorescence intensity to quantify the temporaldynamics of the PTF OCT4 over a cell lifetime. Wecharacterise the internal self-regulation of OCT4 using theHurst exponent and autocorrelation analysis, quantify the intra-cellular fluctuations and consider the diffusive nature of OCT4evolution for individual cells and pairs of their descendants.After cell divisions, OCT4 abundance in the daughter cellsfluctuates sub-diffusively, showing anti-persistent self-regulationwith a Hurst exponent of 0.38. Auto-correlation analysis showsanti-persistence for five hours or longer in 86% of cells, on averagebetween three and 12 hours into the cell cycle. The OCT4fluctuations in five minute intervals follow a Laplace distribution,with BMP4 addition provoking smaller changes and tighter self-regulation, particularly in the differentiated fate group. Thisquantitative framework provides a basis for comparison to otherexperiments, and the development of mathematical models ofpluripotency.
Keywords : human embryonic stem cells, OCT4dynamics, OCT4 regulation, pluripotency
Submitted to:
Phys. Biol.
1. Introduction
Human embryonic stem cells (hESCs) form coloniesthrough repeated mitosis and have the ability todifferentiate into all somatic cell types in the humanbody: the pluripotency property. The pluripotency ofhESCs is their defining characteristic, central to theirtouted applications in drug discovery, regenerative andpersonalised medicine [1, 2, 3, 4, 5, 6]. These promisingclinical applications of hESCs require tight control overcolony pluripotency, homogeneity and differentiationtrajectories in-vitro [7], yet this remains challenging.The control and optimisation of pluripotencyacross colonies is difficult due to the complex inter-regulatory dynamics of pluripotency. At the single-cell level, pluripotency is inherently stochastic. It issuggested that pluripotency is not well defined at thesingle-cell level but is instead a statistical property of acell population [8, 9]. Cell pluripotency is also affectedby many factors: the local environment [10, 11],interactions with neighbours [12, 13], the cell cycle[14] and the substrate [15]. On the colony scale,complex collective effects of pluripotency emerge. Inthe presence of restrictive geometries, differentiatedcells form bands occurring around colony edges [16, 13].Pluripotency maintenance relies on the inter-regulation of pluripotency transcription factors (PTFs):the genes OCT4, SOX2 and NANOG [17, 18, 19]. Af-ter several divisions, PTF fluctuations lead to the es-tablishment of sub-populations with varying pluripo-tency [17]. The differentiation of a stem cell into aspecialised cell is the departure from the pluripotentstate led by PTF destabilisation and their interactionwith chemical signalling pathways [17, 20, 21]. Thisdecision of a stem cell to either remain pluripotent orto differentiate is known as its fate decision. It is un-known how much cell fate decisions are led by inher-ited factors, as opposed to environmental factors andintra-cellular signalling as even clonal (genetically iden-tical) cells under the same conditions make differentfate decisions [22]. Colonies exhibit heterogeneous sub-populations of cells with differing levels of PTF expres-sion [17, 20, 23] which suggests a play-off between thepotentially disruptive single-cell and regulatory com-munity effects [8, 9, 16, 13].On the intra-cellular level, a narrow range of PTFabundance is necessary for maintained pluripotency[24, 25] and even weak fluctuations can bias cellfate decisions in the G1 phase of the cell cycle [26].Furthermore, the PTFs are inherited asymmetricallyas a cell divides, biasing the fate of the daughter cellsand contributing to colony heterogeneity [27, 28, 29]with the decision to differentiate largely determinedbefore any differentiation stimulus is introduced [27].As PTF fluctuations are inherently stochastic [9, 8, 30],it is important to quantify their temporal dynamics a r X i v : . [ q - b i o . CB ] M a y egulatory dynamics of OCT4 in hESCs in-vitro [32] and are informed by recent studies offluctuations of PTFs throughout colonies [27, 26, 9] andthe spatial patterning of differentiation [13, 16]. Manymodels use coupled stochastic differential equations todescribe PTF fluctuations [33, 34, 35] while others usea gene network analysis framework [36, 37] or take amechanistic approach [38].Although the dynamics of OCT4 are complex,affected by many genetic factors and closely regulatedby the other PTFs [17, 39, 21], here we aim to isolateautonomous properties of OCT4 for pluripotent anddifferentiated cells to facilitate the development ofdescriptive mathematical models. This quantificationof OCT4 will provide a basis for identifying systematicsimilarities and differences between PTFs in futureexperiments, and highlights some key indicators ofcell fate. As our quantitative understanding ofPTF regulation increases, more complex regulatoryproperties can be considered to build fundamentalmodels.Here we use the experimental data from Ref.[27] of time-lapse fluorescent measurements of theOCT4-mCherry reporter levels in cells in a growinghESC colony to quantify the dynamics of intra-cellularOCT4. In addition to the OCT4 splitting dynamicsand fluctuations described in Ref. [27], we describequantitatively the fluctuations in OCT4 in relationto cell fate and the addition of the differentiationagent BMP4. We quantify the self-regulation of OCT4through anti-persistence and characterise it within thediffusion framework. The companion manuscript Ref.[40] gives a deeper insight into the spatial distributionof OCT4 in the colony and cell fate transitions, basedon the same experiment. Our quantitative analysisin this manuscript, along with Ref. [27] and [40],provides the basis for developments in mathematicaland statistical models of pluripotency.
2. Materials and Methods
The experiment was carried out by Purvis Lab(University of North Carolina, School of Medicine),and published in Ref. [27]. The OCT4 levels (meanOCT4-mCherry fluorescence intensity) in a human embryonic stem cell (H9) colony were determined andcells were live-imaged for approximately 70 hours. At40 hours the differentiation agent BMP4 was added tothe cells. The cell IDs, ancestries and positions wereextracted along with their OCT4 immuno-fluorescenceintensity values (reported in arbitrary fluorescenceunits a.f.u.). Each cell was classified according to itsfinal fate status as either pluripotent, differentiatedor unknown using expression levels of CDX2. Fullexperimental details are given in Ref. [27]. Thisilluminating study by S. C. Wolff et. al. , provides richopportunities for further quantitative analysis.
The colony begins from 30 cells and grows over 68 hours(817 time frames) to 381 cells, with 1274 cell cycleselapsing within this time. A differentiation agentBMP4 was added to the cells at 40 hours. In Ref. [27]the cells are categorised according to their final cellfate status as pluripotent, differentiated, or unknown(could not confidently be assigned as pluripotent ordifferentiated). The number of cells, N , considered ineach cell fate category, pre- and post-BMP4 is shownin Table 1. N Pre-BMP4 Post-BMP4 All timesPluripotent 96 422 518Differentiated 22 111 133Unknown 112 511 623All fates 230 1044 1274
Table 1: The number of cells, N , in each of the cell fateand pre- and post-BMP4 categories. A post-BMP4 isany cell present at 40 hours or later. There are 1274cells in total.Snapshots of the colony at times T = 0 h, T =20 h, T = 40 h (the time BMP4 is added) and T = 68 h (the final recorded time) colour coded byOCT4 intensity are shown in Figure 1. There isclear spatial patterning of cell fates within the colony,with clustering of pluripotent cells in the centre anddifferentiated cells around the top edge of the colony.Spatial analysis shows that this patterning beginsemerging at around T = 20 hours (20 hours beforeBMP4 addition) for differentiated cells, and at around T = 50 hours (10 hours post BMP4 addition) forpluripotent cells [40]. Although here we focus onquantifying the temporal regulation in OCT4, we mustkeep in mind that there is a spatial correlation betweenthe cell fates.For every cell in the colony there is a correspond-ing time series of the abundance of OCT4 withinthe cell during its lifetime: OCT4( t ), OCT4( t ), ..., egulatory dynamics of OCT4 in hESCs t n ), where t = 0 and t n are the start and end ofthe cell cycle for the cell, respectively. The values of t n range from 15 minutes to 30 hours across the popula-tion. We will use the notation t i to describe time-stepsin terms of the cell cycle and T i for experimental time(between 0 and 68 hours). The OCT4 time series for acell at the beginning of the experiment and its descen-dants are shown in Figure 1(e).An analysis of the number of cells in the colonyover time, N ( T ), is given in the SupplementaryInformation (Figure S1). The whole colony followsexponential growth, with a doubling time of 16 ± . ± . ± . ± .
01 hours for pluripotent, differentiatedand unknown cells respectively. As expected, thepluripotent cells proliferate significantly faster than thedifferentiated cells.
The original dataset is available from Ref. [27],providing cell IDs, cell ancestries, cell positions,cell fates and mean OCT4-mCherry fluorescenceintensities. Cell IDs in this manuscript are consistentwith those in the original dataset.
The quantitative analysis in both Ref. [27] and thismanuscript were performed using MATLAB.
Averaging and errors
When average values are given the type of averaging(mean or median) is specified. For means the errorsare given in the form ± standard deviation (standarderror in the mean). For medians the uncertaintiesquoted represent the lower and upper quartiles or theinterquartile range as specified. Correlation coefficient
The correlation between two OCT4 time series is calcu-lated using Pearson’s correlation coefficient. Note thatthe same conclusion is reached when the Kendall rankcorrelation coefficient is used.
De-trending
We remove trends from the data when it is necessaryto analyse fluctuations about any present trend. Weused MATLAB’s inbuilt function detrend which sub-tracts the best-fit line from the data.
Line fittings
Lines of best fit throughout were calculated using aleast-squares method and the errors given represent the
Time (hours) O C T (e) Figure 1: Snapshots of the colony at (a) T = 0 h,(b) T = 20 h, (c) T = 40 h (at the addition of BMP4)and (d) T = 68 h (final time). The cells are colouredaccording to their OCT4 intensity levels. Note thatthe circles are not indicative of cell or nucleus size. (e)Example OCT4 time series from a cell at the start ofthe experiment (cell ID 25, blue) and its descendentcells.95% confidence interval of the parameters. Statistical testing
To test the null hypothesis that a distribution is Nor-mal, we use both the one-sample Kolmogorov-Smirnovand Shapiro-Wilk tests. To test the null hypothe-sis that two non-parametric distributions are from thesame distribution we use the two-sample Kolmogorov-Smirnov test.
The Laplace distribution
We consider the Laplace distribution, sometimesreferred to as the double exponential distribution, egulatory dynamics of OCT4 in hESCs µ † , b ) to distinguish fromthe usual parameter µ in Normal( µ, σ ). Theparameters can be estimated using the maximumlikelihood estimators ˆ µ † and ˆ b , where ˆ µ † is the samplemedian and ˆ b is the mean absolute deviation from themedianˆ b = 1 N N (cid:88) i − | x i − ˆ µ †| . We use this method of parameter estimation to findLaplace distributions to describe the change in OCT4between time-steps.
Poincar´e maps
Poincar´e maps offer an efficient diagnostic for analysingthe self-similarity of a series. A timed signal is plottedagainst itself after a time delay (here the time delay is5 min) and its scatter pattern reflects the randomnessand variability of the dynamics, thus giving a represen-tation of the correlation between consecutive values ofthe time series. In the Poincar´e plots we colour the val-ues according to their normalised frequency. The data(changes in OCT4) is binned into 100 groups beforethe relative frequency of that group (a value between 0and 1) is calculated. The scatter plot can be quantifiedby fitting an ellipse around the data points and mea-suring the disperson along the major SD1 and minorSD2 axes. Further information on Poincar´e analysis isgiven in Refs. [43] and [44].
The Hurst exponent
The Hurst exponent, H , is a measure of the memory,or the scale of the self-similarity properties of a timeseries. It is defined as E (cid:20) R ( n ) S ( n ) (cid:21) = Cn H as n → ∞ , with R ( n ) the range of the first n cumulative deviationsfrom the mean and S ( n ) their standard deviation. E [ · ]denotes the statistical expectation, n is the number ofdata points in the time series and C is a constant.The quantity R/S is known as the rescaled rangeand measures how the apparent variability changeswith the length of time considered. For a time series X , X , ..., X n , with mean m = n (cid:80) ni =1 X i , R t can becalculated as R i = max( Z , Z , ..., Z t ) − min( Z , Z , ..., Z t ) , where Z t = t (cid:88) i =1 X i − m for t = 1 , , ..., n. Further details on the Hurst exponent, other methodsof estimation and its relation to fractional Brownianmotion can be found in Refs. [45, 46, 47, 48].
Autocorrelation analysis
Autocorrelations were calculated using MATLAB’s autocorr function (Econometrics Toolbox). Theautocorrelation C i of a time series between x t and x t + i for time-lag i is given by C i = 1 T σ T − i (cid:88) t =1 ( x t − x )( x t + k − x ) , where σ is the sample variance of the time series. Thecorrelation time is defined as τ = (cid:82) ∞−∞ C ( t ) dt . Theauto-correlations that occur in the data analysed arewell be described by the function C = cos(2 πt/a ) e − t/b with certain constants a and b [49]. Random walk theory
We apply the theory of random walks to the OCT4time series to test if the OCT4 intensity drifts diffu-sively. Traditionally used for the migration of parti-cles, the mean square displacement (MSD, mean squaredifference or mean square fluctuation) is calculated asMSD = (cid:104)| x ( t ) − x (0) | (cid:105) , where x ( t ) and x (0) are thecurrent (at time t ) and starting positions, and theangular brackets denote the average over all cells in-volved. If the evolution is diffusive (Brownian) thenthe MSD increases linearly with time in the mannerMSD = 2 Dt with a certain constant D known asthe diffusion coefficient, or diffusivity. Sub-diffusionis shown by MSD ∝ t α with α < α >
1. Here we consider the one-dimensional ver-sion, the mean square difference. Further informationon the use of random walks in mathematical biology isgiven in Refs. [50] and [51]. egulatory dynamics of OCT4 in hESCs
3. Results
Firstly, we consider the correlation and similarities inOCT4 expression between pairs of descendant (‘sister’)cells, and proceed to characterise how the drift inOCT4 similarity between these cells emerges overcell lifetimes. We analyse the inherent fluctuationsin OCT4 expression and quantify its self-regulatoryproperties using the Hurst exponent, autocorrelationand diffusion analysis. Where appropriate, we considerthe pluripotent and differentiated cells, pre- and post-BMP4 groups of cells separately to identify anydiagnostic factors.
The analysis in Ref. [27] shows that upon cell divisionthe ratio between the OCT4 values of sister cellsis centred around a 1:1 distribution, meaning thatalthough asymmetric pluripotency inheritance is seen(for example, 38% of divisions occur in the ratio 5:6 ormore extreme), on average sister cells start with similarlevels of OCT4. It is also shown that OCT4 levels aremore similar in closely related cells, i.e., sister cells andcousins cells show significant similarity when comparedwith random pairs of cells.We can consider the strength of the correlation inthe time series for the OCT4 abundance in sister cellsover their whole lifetimes by calculating the correlationcoefficient, ρ . Before calculating the correlation,each OCT4 time series was de-trended to account forany confounding similarities in sister cells that maybe present due to their shared environment. Thedistribution of ρ for pluripotent and differentiated,and pre- and post-BMP4 sister cells are shown inFigure 2(a) and 2(b). Here it is necessary to poolpre- and post-BMP4 cells together for a cell fatecomparison, and vice verse, to keep good statistics.The mean correlations, ρ , are given in Table 2 and showmoderate positive correlations across all categories. ρ Pre-BMP4 Post-BMP4 All timesPluri. - - 0 . ± . . . ± . . . ± . .
03) 0 . ± . .
01) 0 . ± . . Table 2: The mean correlation, ρ ± the standarddeviation (standard error) between pairs of sister cells.There is no difference between cells of differentfates, both with ρ = 0 . ± . ± . ρ for pluripotent and differentiatedcells, shown in Figure 2(a), are the same. Sister cells pre-BMP4 show a weaker correlation than thosepost-BMP4, with ρ = 0 . ± . ρ = 0 . ± . ρ = 0 .
99 and the trend line OCT = (1 ± . .(Note that the labelling of cell 1 and cell 2 is entirelyarbitrary.) By the end of their respective lifetimes, thedistribution spreads, with a correlation of ρ = 0 .
78 anda line of best fit OCT = (0 . ± . .In the next section we will consider the behaviourof OCT4 from the initial point of possible asymmetricinheritance to the final time before mitosis at the endof the cell lifetime and characterise how this drift ofsimilarity in sister cells occurs. In this section we quantify the temporal behaviour ofOCT4 dynamics on the cellular level over the courseof a cell lifetime. We consider the variability betweendiscrete time-steps and quantify the self-regulatorybehaviour of OCT4 using several methods.
To get an overall viewof the OCT4 expression levels, the distributionsof all measured OCT4 values for pluripotent anddifferentiated cells, pre- and post-BMP4 can beconsidered, shown in Figure 3. Pre-BMP4, thedifferentiated cells show a skewed distribution, with anincreased preference towards lower OCT4 expressionsthan the pluripotent cells. This is fitting with the factthat the analysis in Ref. [27] suggests the tendencyto differentiate is largely pre-determined before theaddition of the differentiation stimulus. Post-BMP4,both pluripotent and differentiated cells also show areduction in their OCT4 expression, with the effectseen more strongly in the differentiated cells. It isexpected that the BMP4 causes a reduction in OCT4in the differentiated cells, but it is interesting that thesame effect (to a lesser extent) is also present in thecells which remain pluripotent. The average (median)OCT4 expression levels, OCT4 are shown in Table 3. egulatory dynamics of OCT4 in hESCs -1 0 100.511.5 pd f (a) -1 0 100.511.5 pd f (b) Figure 2: The correlation, ρ , between temporalOCT4 in sister cells where both sisters cells were(a) pluripotent (red unfilled circles) and differentiated(green diamonds) and (b) cells pre-BMP4 (blue filledcircles) and post-BMP4 (orange unfilled squares).OCT4 values for all sister pairs (c) at the start and(d) end of their cell cycles. The lines of best fit(orange solid lines) with standard errors in predictinga future observation (dashed lines) are (c) OCT =(1 ± . with R = 0 .
98 and (d) OCT =(0 . ± . with R = 0 . OCT4 Pre-BMP4 Post-BMP4 All timesPluri. 1500 [1280 1730] 1090 [930 1260] 1160 [980 1380]Diff. 1100 [960 1290] 720 [450 990] 840 [510 1070]All fates 1420 [1190 1670] 1050 [850 1230] 1110 [900 1320]
Table 3: The median OCT4 values, with lower andupper quartiles given as OCT4 [lower quartile upperquartile].
Even smallfluctuations in PTF abundance impact cell fate [26]with both high and low PTF values resulting indifferentiation [24, 25]. Mathematical quantificationof PTF fluctuation will facilitate the description ofpluripotency over discrete time-steps, fitting for time-lapse experiments such as the one considered here[27]. First, we consider the change in the intra-cellular OCT4 abundance between the five minutetime intervals, t , t , ...t n , as ∆OCT4=OCT4( t i )-OCT4( t i − ). Our choice of the five minute intervalsallows for good statistics with 90% of cells having (a) OCT4 pd f -3 (b) OCT4 pd f -3 (c) OCT4 pd f -3 (d) OCT4 pd f -3 Figure 3: The distributions of OCT4 for pluripotentcells (red) (a) pre-BMP4 and (b) post-BMP4, anddifferentiated cells (green) (c) pre-BMP4 and (d) post-BMP4.greater than 50 data points. It is likely that a largeproportion of these individual fluctuations will be dueto experimental noise, but considering all of thesevalues together reveals the average behaviour.The distributions of ∆OCT4 for pluripotent anddifferentiated cells, pre- and post-BMP4 are shownin Figure 4. For both fates across all times, thechange in OCT4 is centred around zero (although theindividual values range from -1300 to 1200). Thismeans that, on average, the change in OCT4 isbi-directional for both pluripotent and differentiatedcells. There is no preference for the abundance toincrease or decrease in a time-step, the fluctuations aresymmetric overall. Interestingly, although symmetric,the distributions are not Normal (confirmed by theKolmogorov-Smirnov and Shapiro-Wilk tests at the95% confidence level) due to a narrower and steeperpeak, shown in Figure 4. A Laplace distribution,Laplace( µ † , b ), better fits the experimental data in allcases, with the parameters given in Table 4. Laplace( µ † , b ) Pre-BMP4 Post-BMP4 All timesPluri. ( − .
7, 52 .
1) ( − .
3, 34 .
7) ( − .
3, 38 . − .
7, 45 .
6) ( − .
6, 23 .
3) ( − .
4, 28 . − .
8, 50 .
8) ( − .
9, 32 .
4) ( − .
4, 28 . Table 4: The parameters from the Laplace( µ † , b )fittings to the ∆OCT4 distributions shown in Figure 4. egulatory dynamics of OCT4 in hESCs (c) -200 0 200 OCT4 pd f (d) -200 0 200 OCT4 pd f (c) -200 0 200 OCT4 pd f (d) -200 0 200 OCT4 pd f Figure 4: Distributions of the change in OCT4between the five minute time frames (∆OCT4) forpluripotent cells (a) pre-BMP4 and (b) post-BMP4,and differentiated cells (c) pre-BMP4 and (d) post-BMP4. Solid lines show the Laplace distributionfittings, Laplace( µ † , b ), with the parameters (a) µ † = − . b = 52 .
1, (b) µ † = − . b = 34 . µ † = − . b = 45 . µ † = − . b = 23 .
3. Dashed lines show the Normal distributionfittings.Post-BMP4 addition, the distributions for bothcell fates become significantly narrower, with the pa-rameter b showing a reduction of 49% for differentiatedcells, and 33% for the pluripotent cells. There is alsoa subtle skew in the differentiated cells towards nega-tive values of ∆OCT4 which is consistent with the factthat the OCT4 levels overall decrease after the BMP4addition. The narrowing of the distributions is dueto a preference of smaller changes in OCT4 in all cellfates provoked by the differentiation agent. This couldbe driven by induced selectivity caused by the BMP4addition (i.e., the BMP4 causes a systematic change,producing a preference for smaller ∆OCT4 values), orit could suggest some collective self-regulation [8]. Fur-ther experiments are needed to investigate if this isa collective behaviour effect, considering the effect ofcolony size. It is expected, since the differentiated cellsare most affected by the BMP4, that this group wouldshow the biggest reduction in variation and thereforethe strongest regulation in their OCT4 values.We can also consider the self-similarity of the Figure 5: Poincar´e maps for the OCT4 signal forpluripotent cells (a) pre-BMP4 and (b) post-BMP4,and differentiated cells (c) pre-BMP4 and (d) post-BMP4. The colour bar shows the normalised relativefrequency of the points.OCT4 series using the Poincar´e map [44, 43]. Foreach cell, its OCT4 time series can be plotted againstitself with one time-step delay, i.e., OCT4( t i ) againstOCT4( t i +1 ) as shown in Figure 5.By assessing qualitatively the shape formed by thereturn map, we observe changes in the distributionof points between pluripotent and differentiated cells,pre- and post-BMP4. Even pre-BMP4 addition, thedifferentiated cells show less variation compared tothe pluripotent cells, with the addition of BMP4exacerbating this effect. Quantitatively these resultscan be described by fitting an ellipse to the shapeformed by the data plots and measuring the dispersionalong the major SD1 and minor SD2 axes, given inTable 5.SD1, SD2 Pre-BMP4 Post-BMP4 All timesPluri. 64, 1430 51, 1150 54, 1320Diff. 54, 1070 33, 990 38, 1090All fates 62, 1440 48, 1170 51, 1330Table 5: Quantitative results for the Poincar´e mapellipse fittings. The major axis (SD1) and minoraxis (SD2) from fitting ellipses to the plots shown inFigure 5. egulatory dynamics of OCT4 in hESCs To investigate the self-regulation and internal memory of OCT4 during a cellcycle, we consider three related approaches, the Hurstexponent, the autocorrelation function and diffusionanalysis.
The Hurst exponent
The Hurst exponent, 0 < H < H = 0 .
5, then the fluctuationsare mutually independent, with the variable just aslikely to increase as decrease at each time-step. Ifthe series is persistent,
H > .
5, then at each time-step the series is more likely to fluctuate in the samedirection as the previous step, i.e., if in the last time-step there was an increase, it is more likely there willbe another increase during the next time-step. Foranti-persistence,
H < .
5, the series is less likely tofluctuate in the same direction as the previous step.The Hurst exponent was calculated for all cellswhich live longer than 50 time frames (4.16 hours).The distributions of all H values for pluripotent anddifferentiated cells, pre- and post-BMP4 are shownin Figure 6. The average Hurst exponents, H ,are given in Table 6 for each group. In all cases,the Hurst exponents are significantly less than 0.5,showing moderate anti-persistence. This shows theself-regulation of OCT4 on the intra-cellular scale: ifthe OCT4 value has just increased, it is more likely tonext decrease, and vice versa. This is the case acrosseach cell fate group. Although the means are withinerrors of one another, the Kolmogorov-Smirnov testsreject the null hypothesis that the distributions of H for pluripotent and differentiated pre-BMP4 cells arethe same at the 95% level. There is no significantdifference in H before and after the BMP4 additionfor both cell fates (confirmed by the Kolmogorov-Smirnov test at the 95% level) suggesting this aspectof the self-regulatory behaviour is inherent within thecells and unchanged by the differentiation stimulus.This quantification via the Hurst exponent is directlytransferable to use in fractional Brownian motionmodelling methods [45, 46, 47, 48]. (a) H pd f (b) H pd f (c) H pd f (d) H pd f Figure 6: The Hurst exponent, H for pluripotentcells (red) (a) pre-BMP4 and (b) post-BMP4 anddifferentiated cells (green) (c) pre-BMP4 and (d) post-BMP4. The black lines show H = 0 .
5, the value forBrownian fluctuations. H Pre-BMP4 Post-BMP4 All timesPluri. 0.37 (0.08 0.008) 0.37 (0.09 0.004) 0.37 (0.09 0.004)Diff. 0.42 (0.08 0.02) 0.39 (0.09 0.009) 0.40 (0.09 0.008)All fates 0.38 (0.08 0.007) 0.38 (0.09 0.004) 0.38 (0.09 0.004)
Table 6: The mean Hurst exponent H with (standarddeviation, standard error) for all cell categories. Autocorrelation
The anti-persistence can be furtherexplored by considering the autocorrelation of theOCT4 time series. The autocorrelation is thecorrelation of a time series with itself at increasingtime lags, hence − ≤ C ≤ C = 0signifies no correlation, C <
C > egulatory dynamics of OCT4 in hESCs C (a) Time lag (hours)
Time (hours) O C T (d) C (b) Time lag (hours)
Time (hours) O C T (e) C (c) Time lag (hours)
Time (hours) O C T (f) Figure 7: Typical autocorrelations showing (a) aperiod of anti-correlation before settling at zerocorrelation (seen in 51% of cells), (b) two periods ofanti-correlation followed by correlation (seen in 28% ofcells) and (c) a period of anti-correlation followed by aperiod of correlation (seen in 14% of cells). The panels(d)-(f) show the OCT4 variation in time for these cellsrespectively. The average behaviour is similar to thatin (a).the autocorrelation settles at zero. There are, however,other behaviours evident. Some cells show several lagintervals of anti-correlation, as in Figure 7(b) (Cell ID14), with others showing a positive correlation at alonger time lag before settling at zero, as in Figure 7(c)(Cell ID 43). The corresponding time series of OCT4for each example cell are shown in Figure 7(d-f). Theintervals of positive correlation are visible as trends inOCT4 (either continued increase or decline), with anti-correlation visible as fluctuations about a horizontalline. Anti-correlation for a time lag of at least one hourduration is seen in 99% (1255/1274) of cells, and forat least five hours in 86% (1090/1274) of cells. Ofthe cells with at least one hour anti-correlation visible,44% show a second period of correlation near the endof their lifetimes (as in Figure 7(a) and (b)). Out ofthese, 65% show one period of anti-correlation (as inFigure 7(c)), 31% two periods (as in Figure 7(b)), andthe remaining 4% three or more. For the 57% of cellswith no second period of correlation, 90% of cells showone period, 8% two periods and 2% three or moreperiods of anti-correlation. There is no correlation between the number of periods of anti-correlation andcell fate or the cell’s average position in the colony.The first time anti-correlation occurs, t AC , can becalculated for each individual cell. The distribution of t AC for cells with at least one hour anti-correlation isshown in the Supplementary Information, Figure S2and reveals the critical cell cycle time in which it firstoccurs. In all cells with anti-correlation, it has begunby 8 hours into the cell cycle (just over half a cellcycle [27]), suggesting that before they reach the latterhalves of their lifetimes the internal self-regulation ofOCT4 begins. This could be due to the memory effectsor the down-regulation of the PTF which occurs priorto mitosis [52, 53].The periodic nature and decay of the auto-correlation can be captured by the function C =cos(2 πt/a ) e − t/b [49] (note that this periodicity in theautocorrelation does not necessarily imply periodicityin the time series). These fittings are shown in the Sup-plementary Information, Figure S3, for 25 randomlyselected cells in the colony. This quantifies the tem-poral, periodic decay in the autocorrelation, with theparameter a representing the time-scale of the period-icity, and b the time-scale of the decay (the correlationdecay time). Histograms of a and b for all 1274 cellsare shown in Figure 8. Both distributions are skewed,with medians of 11.7 h and 3.0 h, and 90th percentilesof 30 h and 7 h for a and b respectively. This quanti-fies the characteristic time-scale of the periodicity andthe correlation decay time as less than 7 hours in 90%of cases. For consistency, the parameters split by cellfate pre- and post-BMP4 are given in the Supplemen-tary Information, Figure S4 and S5.The correlation time is defined as τ = (cid:82) ∞−∞ C ( t ) dt ,with a mean correlation time across all cells of τ ≈ ± .
002 h. The distribution of all correlation times isshown in the Supplementary Information, Figure S2.We can identify the average behaviour byconsidering all autocorrelations for all cells. Themean (and standard deviation) and median (andinterquartile range) autocorrelations C for all cellsis shown in Figure 8(c). Notably the mean andmedian are comfortably within errors of each other andthe autocorrelation is robust to the chosen averagingmethod. The average autocorrelation decreases to zeroat around three hours, followed by a period of negativeautocorrelations indicative of anti-persistent behaviourbetween approximately three and 12 hours. By 13hours, the average autocorrelation settles at zero,showing no internal memory past this time. Theseobservations are robust to cell fate and the equivalentautocorrelations for pluripotent and differentiated cellsare shown in the Supplementary Information, FigureS6. This shows that during a cell cycle, there is long-term memory in the OCT4 expression up to around egulatory dynamics of OCT4 in hESCs (a) a (hours) pd f (b) b (hours) pd f Time (hours) -0.500.51 C (c) Figure 8: The distributions of the parameters (a) a and (b) b from all C = cos(2 πt/a ) e − t/b autocorrelationfittings. The parameter estimates using the mean andmedian autocorrelations for all cells are a = 11 . ± . . ± .
69 and b = 2 . ± .
23 and 2 . ± . πt/a ) e − t/b , as the full scaleof the anti-correlation is not captured, shown in theSupplementary Information, Figure S6. Diffusion analysis
A further method of quantifyingthe internal regulation of OCT4 is to consider thediffusive behaviour of the time series. The theory ofdiffusion and random walks is widely used across manybiological applications, including stem cells and so it isimportant to quantify the OCT4 behaviour within thisframework [50, 51, 54, 55, 56, 57].After the cell division with asymmetric inheritanceof OCT4 there is a short period of increasedfluctuations [27]. Here we consider each OCT4 timeseries from half an hour after cell division to explore theproperties of these fluctuations. Each cell has an OCT4value at the start of its lifetime, denoted OCT . Themean square difference of OCT4 over time, MSD( t ), can be calculated as (cid:104)| OCT4( t ) − OCT4 | (cid:105) , where theangular brackets denote the average across all cells inthe group considered. The MSD for pluripotent anddifferentiated cells, pre- and post-BMP4 between 0 and12 hours is shown in Figure 9. For pluripotent cells,both pre- and post-BMP4, the distinct sub-diffusivebehaviour of the MSD is visible, with MSD = βt α , α <
1. The parameters α and β are shown in Table 7.The differentiated cells do not follow a power lawrelationship, pre- or post-BMP4, but the limiting ofthe MSD can still be seen from around 2 hours pre-BMP4 and from 9 hours post-BMP4.Pre-BMP4 Post-BMP4 α . ± .
03 0 . ± . β ± ± βt α fittings forpluripotent cells.This sub-diffusivity is consistent with the anti-persistence revealed using the Hurst exponent andautocorrelation. In agreement with those results, onaverage, the intra-cellular OCT4 abundance behaves ina sub-diffusive manner throughout a cell lifetime. Thishas a knock-on effect for the relationship between sistercell OCT4 which is presented in the SupplementaryInformation (Figure S7-9).This further quantifies the self-regulatory be-haviour of OCT4 within the diffusion framework, a fun-damental starting point for many mathematical mod-els. The anti-persistence of OCT4 suggests possibilitiesfor mathematical modelling methods to capture the in-ternal regulation of pluripotency, including fractionalBrownian motion and correlated random walk theory. We have quantified the behaviour of OCT4 using avariety of mathematical techniques, some of whichcould provide a non-invasive diagnostic tool to identifypluripotent and differentiated cells. To illustrate this,we will use the unknown cells (unable to be classifiedas either pluripotent or differentiated) and comparetheir time series parameters to the pluripotent anddifferentiated cells.Firstly, the distribution of all OCT4 values for theunknown cells lies between that for the pluripotentand differentiated cells, shown in Figure 10(a).This is unsurprising as these cells had an OCT4expression (along with a CDX2 expression) that didnot correspond to either cell fate. Having measuredthe OCT4 time series for enough time steps to get adistribution of OCT4, comparison could be made toidentify whether it bests corresponds to that of the egulatory dynamics of OCT4 in hESCs Time (hours) M S D (a) Time (hours) M S D (b) Figure 9: The MSD for pluripotent (red solid line)and differentiated (green dashed line) between 0 and12 hours (a) pre-BMP4 and (b) post-BMP4. Theblack dashed lines show the fits MSD = βt α with theparameters (a) α = 0 . ± . β = 42000 ± α = 0 . ± . β = 35000 ± µ † = − . b =27, both between their pluripotent and differentiatedcounterparts. The distributions pre-BMP4 aren’tdifferent enough to distinguish between the fates,shown in the Supplementary Information Figure S10.There is no distinguishable difference betweenthe cell fates using an autocorrelation analysis or theHurst exponent (with H within errors for all cell fates,including the unknown cells). However, the MSD plotsto identify sub-diffusion show a significant differencebetween the fates, with the MSD for pluripotent cells (a) OCT4 pd f -3 (b) -200 0 200 OCT4 pd f Time (hours) M S D (c) Figure 10: (a) The distribution of all OCT4 values forunknown cells pre-BMP4. (b) The distribution of allchanges in OCT4, ∆OCT4, for all unknown cells post-BMP4. In both cases the corresponding distributionsfor pluripotent and differentiated cells as a red solidline and green dashed line, respectively. (c) The MSDfor unknown cells pre-BMP4 (blue solid line) and post-BMP4 (orange dashed line) with standard error errorbars.well described by a power law relationship. The MSDfor the unknown cells in shown in Figure 10(c). Pre-BMP4, the MSD shows significant sub-diffusivity in thefirst 4 hours, but can not be described by a power law,unlike for pluripotent cells. The MSD also shows moresimilarities with the differentiated cells post-BMP4,with a more linear MSD which levels off at around 11hours.These results show that overall, the unknowncells behaviour lies between that of the pluripotentand differentiated cell fates. The distinct differencesbetween cell fates, seen in the OCT4 distributions,the change in OCT4 post-BMP4 and the sub-diffusive nature of the MSD could provide non-invasivediagnostic tools to identify cell fates. egulatory dynamics of OCT4 in hESCs
4. Discussion
Promising clinical applications of hESCs require tightcontrol over the pluripotency of hESC colonies. Ithas been shown that even small PTF fluctuations canbias cell fate decisions and that PTFs are inheritedasymmetrically upon cell division [26, 27, 28, 29].It is therefore necessary to quantify the dynamicsof key PTFs to further our understanding of howpluripotency is regulated and assist in the developmentof mathematical modelling. Thorough quantificationalso provides the basis for experimental comparisons,and the identification of systematic and universalbehaviours. Here we have used a published data setfrom Ref. [27] to analyse and quantify the dynamics ofthe pluripotency transcription factor OCT4.The colony considered here grows exponentially,with changing proportions of pluripotent, differenti-ated and unknown cells. Snapshots of the colony showsome spatial patterning of the OCT4 abundance (Fig-ure 1), with higher levels of expression of OCT4 vis-ible in cells clustered in the colony centre. A spatialanalysis of the colony can be found in the companionmanuscript, Ref. [40]. Here we have focused on thequantification of the temporal dynamics of OCT4.Time-lapse experiments such as the one consideredhere [27] provide opportunities for the quantificationof temporal PFT regulation which can be comparedto, and enhance, current biological knowledge. Forexample, a sharp decline in OCT4 levels occurringbefore cell division is noted in Ref. [27], in keepingwith the transcription factor down regulation knownto occur before mitosis [52, 53]. This phenomena canbe quantified, with the decrease in OCT4 beginning,on average, 35 minutes (0.58 hours) before celldivision, lasting for 15 minutes (0.25 hours), andshowing a reduction of 22% from the interphaseOCT4 expression. The OCT4 expression levels recoveras mitosis occurs and the cycle repeats for thecell’s descendants consistent with experimental resultsshowing OCT4 resets on re-entry to the G1 phase[58, 59]. This is shown for all cells before BMP4addition in the Supplementary Information, FigureS11. Ref. [27] reveals that sister cells show moreclosely related OCT4 values than pairs of randomcells. Here we take this a step further by quantifyingtheir temporal dynamics in relation to one another.Taking into account any common trends which mayaffect both cells due to their shared environment, thesister cells before BMP4 show a moderate correlationwith each other with a correlation coefficient of 0.5.This is reduced to a slight correlation for pairs thatexist after the BMP4 addition (0.3). The fact thatthese correlations still occur after de-trending furtherhighlights the inherent similarities between related cells. We then consider the OCT4 behaviour over celllifetimes to explore the manner in which this drift insimilarity occurs. The behaviour is summarised in theschematic in Figure 11.Stochastic fluctuations in OCT4 have been shownto bias cell fate [26] with evidence of asymmetric noiseleading to noise-mediated cell plasticity [30]. Herewe see the change in OCT4 between each 5 minutetime interval is isotropic, with an average of zero. Anatural assumption in model development would beto simulate this symmetric time-step change in OCT4with a Normal distribution, however the distributionof all these changes best fits a Laplace distribution.Further experimental data are needed to confirm this isa robustly appropriate choice, elucidate the parametersfor other experimental conditions and investigate howthis is affected by cell-cell interactions. Note that thisallows us to capture the nature of the variation inOCT4 only and further aspects of the behaviour needto be considered to fully describe the OCT4 regulationover time.Although this shows that overall, positive changesin OCT4 are just as likely to occur as negative ones, itdoes not reveal anything about the temporal nature ofthese fluctuations and hence any correlation propertieswhich may be evident over time (for example, all thepositive changes in OCT4 could come one after theother, followed by all the negative changes, it doesnot mean that a positive change is necessarily followedby a negative change). There is also a difference inthese fluctuations after the differentiation agent, withthe addition of BMP4 provoking tighter self-regulationacross all cell fates. Further experiments are needed toinvestigate whether this self-regulation is a collectivebehaviour effect.A significant finding of this analysis is thequantification of the self-regulatory properties of OCT4within cells. An autocorrelation analysis, alongwith the calculation of the Hurst exponent of 0.38shows significant anti-persistence, in keeping with theregulation of PTFs [17, 21, 37]. Throughout the colonygrowth, anti-correlation of at least five hours is seen in86% of cells (with no significant difference between thecell fates, suggesting OCT4 regulation is quantitativelycomparable in pluripotent and differentiated cells), andon average occurs between 3 and 12 hours into a cell’slifetime. This is further illustrated by consideringthe behaviour of the cells in the diffusion framework,with cells across all fates showing significant sub-diffusivity. For pluripotent cells, the sub-diffusivitycan be described by a power-law relationship. Thesub-diffusivity allows another characterisation of theself-regulation of OCT4 expression and provides aquantitative starting point for the mathematicalmodelling of OCT4 time series. egulatory dynamics of OCT4 in hESCs Inheritance≈1:1
𝐎𝐂𝐓𝟒 𝚫𝐎𝐂𝐓𝟒 ~ Sub-diffusive
Pluri. 1500 L(-0.7, 52.1) MSD=
Diff. 1100 L(-0.7,45.6) -
OCT4( t )OCT4( t ) (a) (c)(d) 𝐎𝐂𝐓𝟒 𝚫𝐎𝐂𝐓𝟒 ~ Sub-diffusive
MSD = 35000𝑡
720 L(-2.6,23.3) -BMP4 (b)
Figure 11: An illustration of the dynamics in OCT4 over a cell lifetime. (a) OCT4 is split, possibly asymmetricallybut on average in a 1:1 ratio [27] before fluctuating in a sub-diffusive manner (with a power law relationship forpluripotent cells), resulting in (b) more variation in sister cells at the end of their lifetimes (over an average of14 hours [27]). This can result in cells of different fates: (c) pluripotent and (d) differentiated.
Implications
The quantification of the temporal fluctuations inPTFs is essential for experimental comparisons and thedevelopment of mathematical models of pluripotency.We have further developed the experimental analysisprovided in Ref. [27] to quantify the self-regulation ofOCT4 over a cell lifetime. These results can be usedto guide mathematical modelling; informing equationchoices, revealing properties of the inherent chemistryand identifying universal behaviours, and facilitateexperimental diagnostic comparisons (as we consideredin Section 3.3 for the unknown cell fate category).The distributions of OCT4 provide a quantitativecomparison for future experiments, characterising thespread and skew of OCT4 expression. Similarly, theLaplace distribution fittings to the change in OCT4between the five minute time intervals can be a directinput into computation models requiring ∆OCT4 asa probabilistic parameter. The distributions canalso be used for model verification. Regardlessof the modelling method used, the distributions ofOCT4 and ∆OCT4 produced via a simulation canbe compared, both qualitatively, and quantitativelyby the parameters of the distribution, to theseexperimental results.The Hurst exponent quantifies the level of anti-persistence in the OCT4 regulation. Broadly, theidentification of a Hurst exponent which is notBrownian suggests the use of specific equations, i.e.,the equation needs to be able to produce this anti-persistence. Furthermore, it can be a direct inputinto some stochastic modelling equations in which theparameter, H , is required, e.g., fractional Brownianmotion [45, 46, 47, 48].The results suggest that the anti-correlation is an inherent property of the cells, across all cell fates.Further experiments are needed to clarify that thisis the case under different experimental conditionsbut this provides a further quantitative statistic forcomparisons. The identification of this systematicproperty has implications for the underlying stochasticchemistry of the OCT4 regulation and can be used toinform chemical models of PTF regulation, often basedon the Hill equations [33, 60].Furthermore, the sub-diffusive nature of thetime series for pluripotent cells allows the completecharacterisation using only two summary parameters, α and β . This highlights the presence of a universalbehaviour in the cells, and can be used as a directparameter input into Brownian (and the fractional andgeometric extensions) computational models allowingthe future behaviour of the OCT4 regulation to bepredicted statistically.The experiment in Ref. [27] has led to a richanalysis, allowing us to establish the language throughwhich to quantitatively compare this experiment toothers and to guide mathematical modelling choices.In general, this highlights the need for further temporalexperimental data on OCT4 and other transcriptionfactors. Author Contributions
Author contributions for the original experiment givenin [27]. L.E.W. and S.O.F. analysed the data andprepared the figures. L.E.W., S.O.F., I.N., M.L.,N.G.P. and A.S. wrote the manuscript. egulatory dynamics of OCT4 in hESCs Acknowledgements
SOF acknowledges financial support from the ConsejoNacional de Ciencia y Tecnolog´ıa (CONACyT, Mex-ico) for the grant CVU-174695. IN acknowledges thegrant from the Russian Government 641 Program forthe recruitment of the leading scientists into 641 Rus-sian Institution of Higher Education 14.w03.31.0029and RFFI project grant number 20-015-00060. ML ac-knowledges BBSRC UK (BB/I020209/1) for providingfinancial support for this work.
References [1] A. D. Ebert and C. N. Svendsen. Human stem cells anddrug screening: opportunities and challenges.
Nat. Rev.Drug Discov. , 9(5):367–372, 2010.[2] Z. Zhu and D. Huangfu. Human pluripotent stem cells: anemerging model in developmental biology.
Development ,140(4):705–717, 2013.[3] Y. Avior, I. Sagi, and N. Benvenisty. Pluripotent stem cellsin disease modelling and drug discovery.
Nat. Rev. Mol.Cell Biol , 17(3):170–182, 2016.[4] D. Ilic and C. Ogilvie. Concise Review: Human embryonicstem cells - what have we done? what are we doing?where are we going?[5] G. Shroff, J. D. Titus, and R. Shroff. A review of theemerging potential therapy for neurological disorders:human embryonic stem cell therapy.
Am. J. Stem Cells ,6(1):1, 2017.[6] A. Trounson and N. D. DeWitt. Pluripotent stem cellsprogressing to the clinic.
Nat. Rev. Mol. Cell Biol ,17(3):194, 2016.[7] C. L. Bauwens, R. Peerani, S. Niebruegge, K. A.Woodhouse, E. Kumacheva, M. Husain, and P. W.Zandstra. Control of human embryonic stem cellcolony and aggregate size heterogeneity influencesdifferentiation trajectories.
Stem Cells , 26(9):2300–2310,2008.[8] B. D. MacArthur and I. R. Lemischka. Statisticalmechanics of pluripotency.
Cell , 154(3):484 – 489, 2013.[9] M-E Torres-Padilla and I. Chambers. Transcription factorheterogeneity in pluripotent stem cells: a stochasticadvantage.
Development , 141(11):2173–2181, 2014.[10] E. Shuzui, M. Kim, and M. Kino-oka. Anomalouscell migration triggers a switch to deviation from theundifferentiated state in colonies of human inducedpluripotent stems on feeder layers.
J. Biosci. Bioeng. ,127(2):246 – 255, 2019.[11] R. Stadhouders, G. J. Filion, and T. Graf. Transcriptionfactors and 3D genome conformation in cell-fatedecisions.
Nature , 569(7756):345–354, 2019.[12] A. Nemashkalo, A. Ruzo, I. Heemskerk, and A. Warmflash.Morphogen and community effects determine cell fatesin response to bmp4 signaling in human embryonic stemcells.
Development , 144(17):3042–3053, 2017.[13] K. A. Rosowski, A. F. Mertz, S. Norcross, E. R. Dufresne,and V. Horsley. Edges of human embryonic stemcell colonies display distinct mechanical properties anddifferentiation potential.
Sci. Rep. , 5:14218, Sep 2015.[14] S. Pauklin and L. Vallier. The cell-cycle state of stem cellsdetermines cell fate propensity.
Cell , 155(1):135 – 147,2013.[15] N. S. Hwang, S. Varghese, and J. Elisseeff. Controlleddifferentiation of stem cells.
Adv. Drug Deliv. Rev. ,60(2):199 – 214, 2008. Emerging Trends in Cell-BasedTherapies. [16] A. Warmflash, B. Sorre, F. Etoc, E. D. Siggia, and A. H.Brivanlou. A method to recapitulate early embryonicspatial patterning in human embryonic stem cells.
Nat.Methods , 11(8):847–854, 2014.[17] M. Li and J. C. Izpisua Belmonte. Deconstructing thepluripotency gene regulatory network.
Nat. Cell Biol. ,20(4):382–392, 2018.[18] L. A. Boyer, T. I. Lee, M. F. Cole, S. E. Johnstone, S. S.Levine, J. P. Zucker, M. G. Guenther, R. M. Kumar,H. L. Murray, R. G. Jenner, D. K. Gifford, D. A. Melton,R. Jaenisch, and R. A. Young. Core transcriptionalregulatory circuitry in human embryonic stem cells.
Cell ,122(6):947—956, September 2005.[19] I. Chambers and S. R. Tomlinson. The transcriptionalfoundation of pluripotency.
Development , 136(14):2311–2322, 2009.[20] R. M. Kumar, P. Cahan, A. K. Shalek, R. Satija,A. DaleyKeyser, H. Li, J. Zhang, K. Pardee, D. Gennert,J. J. Trombetta, T. C. Ferrante, A. Regev, G. Q.Daley, and J. J. Collins. Deconstructing transcriptionalheterogeneity in pluripotent stem cells.
Nature ,516(7529):56–61, 2014.[21] Z. Wang, E. Oron, B. Nelson, S. Razis, and N. Ivanova.Distinct lineage specification roles for NANOG, OCT4,and SOX2 in human embryonic stem cells.
Cell StemCell , 10(4):440 – 454, 2012.[22] O. Symmons and A. Raj. What’s luck got to do withit: single cells, multiple fates, and biological non-determinism.
Mol. Cell , 62(5):788 – 802, 2016.[23] A. M. Singh, J. Chappell, R. Trost, L. Lin, T. Wang,J. Tang, H. Wu, S. Zhao, P. Jin, and S. Dalton. Cell-cycle control of developmentally regulated transcriptionfactors accounts for heterogeneity in human pluripotentcells.
Stem Cell Reports , 1(6):532 – 544, 2013.[24] H. Niwa, J. Miyazaki, and A. G. Smith. Quantitative ex-pression of Oct-3/4 defines differentiation, dedifferentia-tion or self-renewal of ES cells.
Nat. Genet. , 24(4):372–376, 2000.[25] J. L. Kopp, B. D. Ormsbee, M. Desler, and A. Rizzino.Small increases in the level of Sox2 trigger thedifferentiation of mouse embryonic stem cells.
Stem cells ,26(4):903–911, 2008.[26] D. Strebinger, C. Deluz, E. T. Friman, S. Govindan, A. B.Alber, and D. M. Suter. Endogenous fluctuations of oct4and sox2 bias pluripotent cell fate decisions.
Mol. Syst.Biol. , 15(9):e9002, 2019.[27] S. C. Wolff, K. M. Kedziora, R. Dumitru, C. D. Dungee,T. M. Zikry, A. S. Beltran, R. A. Haggerty, J. Cheng,M. A. Redick, and J. E. Purvis. Inheritance of oct4predetermines fate choice in human embryonic stem cells.
Mol. Syst. Biol. , 14(9):e8140, 2018.[28] M. Skamagki, K. B. Wicher, A. J., S. Ganguly, andM. Zernicka-Goetz. Asymmetric localization of Cdx2mRNA during the first cell-fate decision in early mousedevelopment.
Cell Reports , 3(2):442 – 457, 2013.[29] W-W. Tee and D. Reinberg. Chromatin features and theepigenetic regulation of pluripotency states in ESCs.
Development , 141(12):2376–2390, 2014.[30] W. R. Holmes, N. S. Reyes de Mochel, Q. Wang, H. Du,T. Peng, M. Chiang, O. Cinquin, K. Cho, and Q. Nie.Gene expression noise enhances robust organization ofthe early mammalian blastocyst.
PLoS Comput. Biol. ,13(1):1–23, 01 2017.[31] L. E. Wadkin, S. Orozco-Fuentes, I. Neganova, M. Lako,A. Shukurov, and N. G. Parker. The recent advances inthe mathematical modelling of human pluripotent stemcells.
SN Applied Sciences , 2(2):276, 2020.[32] Pınar Pir and Nicolas Le Nov`ere. Mathematical modelsof pluripotent stem cells: at the dawn of predictiveregenerative medicine. In
Systems Medicine , pages 331– egulatory dynamics of OCT4 in hESCs PLoS Comput. Biol. , 2(9):1–13, 092006.[34] I. Glauche, M. Herberg, and I. Roeder. Nanog variabilityand pluripotency regulation of embryonic stem cells -insights from a mathematical model analysis.
PLoS One ,5(6):1–12, 06 2010.[35] M. Herberg and I. Roeder. Computational modellingof embryonic stem-cell fate control.
Development ,142(13):2250–2260, 2015.[36] H. Xu, Y-S. Ang, A. Sevilla, I. R. Lemischka, andA. Ma’ayan. Construction and validation of a regulatorynetwork for pluripotency and self-renewal of mouseembryonic stem cells.
PLoS Comput. Biol. , 10(8):1–14,08 2014.[37] I. R. Akberdin, N. A. Omelyanchuk, S. I. Fadeev, N. E.Leskova, E. A. Oschepkova, F. V. Kazantsev, Y. G.Matushkin, D. A. Afonnikov, and N. A. Kolchanov.Pluripotency gene network dynamics: System views fromparametric analysis.
PLoS One , 13(3):1–24, 03 2018.[38] D. Auddya and B. J. Roth. A mathematical description ofa growing cell colony based on the mechanical bidomainmodel.
J. Phys. D. Appl. Phys. , 50(10):105401, 2017.[39] Yasmin Babaie, Ralf Herwig, Boris Greber, Thore CBrink, Wasco Wruck, Detlef Groth, Hans Lehrach, TomBurdon, and James Adjaye. Analysis of Oct4-dependenttranscriptional networks regulating self-renewal andpluripotency in human embryonic stem cells.
Stem cells ,25(2):500–510, 2007.[40] S. Orozco-Fuentes, L. E. Wadkin, I. Neganova, M. Lako,R. A. Barrio, A. W. Baggaley, A. Shukurov, and N. G.Parker. Spatio-temporal analyses of OCT4 expressionand fate transitions in human embryonic stem cells. bioRxiv , 2020.[41] P. N. Ghule, R. Medina, C. J. Lengner, M. Mandeville,M. Qiao, Z. Dominski, J. B. Lian, J. L. Stein, A. J.van Wijnen, and G. S. Stein. Reprogramming thepluripotent cell cycle: Restoration of an abbreviated G1phase in human induced pluripotent stem (iPS) cells.
J.Cell. Physiol. , 226(5):1149–1156, 2011.[42] L. E. Wadkin, S. Orozco-Fuentes, I. Neganova, S. Bojic,A. Laude, M. Lako, N. G. Parker, and A. Shukurov.Seeding hESCs to achieve optimal colony clonality.
Sci.Rep. , 9:15299, 2019.[43] M. Fishman, F. J. Jacono, S. Park, R. Jamasebi,A. Thungtong, K. A. Loparo, and T. E. Dick. A methodfor analyzing temporal patterns of variability of a timeseries from Poincar´e plots.
J. Appl. Physiol. , 113(2):297–306, 2012.[44] A. Burykin, M. D. Costa, L. Citi, and A. L. Goldberger.Dynamical density delay maps: simple, new method forvisualising the behaviour of complex systems.
BMCMed. Inform. Decis. Mak. , 14(1):6, 2014.[45] B. B. Mandelbrot and J. W. Van Ness. Fractional brownianmotions, fractional noises and applications.
SIAM Rev ,10(4):422–437, 1968.[46] J. Mielniczuk and P. Wojdy(cid:32)l(cid:32)lo. Estimation of Hurstexponent revisited.[47] L. Lacasa, B. Luque, J. Luque, and J. C. Nuno. Thevisibility graph: A new method for estimating the Hurstexponent of fractional Brownian motion. 86(3):30001,2009.[48] J. Barunik and L. Kristoufek. On Hurst exponentestimation under heavy-tailed distributions.
PHYSICAA , 389(18):3844–3855, 2010.[49] A. A. Sveshnikov, I. N. Sneddon, and M. Stark.
AppliedMethods of the Theory of Random Functions . ISSN.Elsevier Science, 1966. [50] J. D. Murray.
Mathematical Biology I. An Introduction ,volume 17 of
Interdisciplinary Applied Mathematics .Springer, New York, 3 edition, 2002.[51] E. A. Codling, M. J. Plank, and S. Benhamou. Randomwalk models in biology. 6(25):813–834, 2008.[52] K. S. Zaret. Genome reactivation after the silence inmitosis: recapitulating mechanisms of development?
Dev. Cell , 29(2):132–134, 2014.[53] N. Festuccia, I. Gonzalez, N. Owens, and P. Navarro.Mitotic bookmarking in development and stem cells.
Development , 144(20):3633–3645, 2017.[54] L. Li, B. H. Wang, S. Wang, L. Moalim-Nour, K. Mohib,D. Lohnes, and L. Wang. Individual cell movement,asymmetric colony expansion, rho-associated kinase, andE-cadherin impact the clonogenicity of human embryonicstem cells.
Biophys. , 98:2442 – 2451, 2010.[55] P. Wu, A. Giri, S. X. Sun, and D. Wirtz. Three-dimensional cell migration does not follow a randomwalk. 111(11):3949–3954, 2014.[56] L. E. Wadkin, L. F. Elliot, I. Neganova, N. G. Parker,V. Chichagova, G. Swan, A. Laude, M. Lako, andA. Shukurov. Dynamics of single human embryonic stemcells and their pairs: a quantitative analysis.
Sci. Rep. ,7(1):1–12, 2017.[57] L. E. Wadkin, S. Orozco-Fuentes, I. Neganova, G. Swan,A. Laude, M. Lako, A. Shukurov, and N. G. Parker.Correlated random walks of human embryonic stem cellsin vitro.[58] Jihoon Shin, Tae Wan Kim, Hyunsoo Kim, Hye Ji Kim,Min Young Suh, Sangho Lee, Han-Teo Lee, SojungKwak, Sang-Eun Lee, Jong-Hyuk Lee, Hyonchol Jang,Eun-Jung Cho, and Hong-Duk Youn. Aurkb/pp1-mediated resetting of oct4 during the cell cycledetermines the identity of embryonic stem cells. eLife ,5:e10877, 2016.[59] Hye Ji Kim, Jihoon Shin, Sangho Lee, Tae Wan Kim,Hyonchol Jang, Min Young Suh, Jae-Hwan Kim, In-Young Hwang, Deog Su Hwang, Eun-Jung Cho, andHong-Duk Youn. Cyclin-dependent kinase 1 activitycoordinates the chromatin associated state of Oct4during cell cycle in embryonic stem cells.
Nucleic AcidsResearch , 46(13):6544–6560, 06 2018.[60] V. Likhoshvai and A. Ratushny. Generalized Hill functionmethod for modeling molecular processes.