Topological transitions, turbulent-like motion and long-time-tails driven by cell division in biological tissues

The cell, an active agent, is the basic unit of life. The motion of a single cell has been studied intensively in the past few decadeslauffenburger1996cell ; horwitz1999cell ; bray2000cell ; ponti2004two ; li2009actin ; raab2012crawling ; novikova2017persistence ; van2018mechanoreciprocity . However, cells interact with each other and move in a collective manner, which plays an essential role in many biological processes, such as embryonic morphogenesis, organ development, wound healing, and spreading of cancerfriedl2009collective ; rorth2009collective ; poujade2007collective . The collective movement of cells using physical models, spanning from sub-cellular to supra-cellular length scaleshakim2017collective ; camley2017physical ; alert2020physical , gives rise to unexpected spatial and temporal correlations.

Since the discovery that flow-induced stretching of a long polymer in a high viscosity solution could produce flow patterns resembling that found in fully developed turbulence groisman2000elastic , a number of studies have shown that turbulent-like features emerge in a variety of biological systems alert2022active ; dombrowski2004self ; rossen2014long ; bratanov2015new ; wensink2012meso . In these examples, the complex flow patterns are thought to rise from active forces, and hence the collective behavior is referred to as active turbulence, even though the scaling behavior is not universal. Hydrodynamic models, that build on the Toner-Tu equations toner1998flocks ; toner2005hydrodynamics ; wensink2012meso ; bratanov2015new and particle based simulations for self-propelled rods wensink2012meso have generated turbulent-like patterns. In addition, numerical solutions and theoretical arguments of active nematic liquid crystal models have been used to show that vortex formation and turbulent-like behavior emerges giomi2015geometry ; thampi2015intrinsic ; alert2020universal . In all these models, active forces, the cause for turbulent-like behavior, is introduced phenomenologically (by “hand” as it were).

Recent studiesrossen2014long ; doostmohammadi2015celebrating suggest that cell division and apoptosis lead to turbulent-like velocity fields with long-range vortex structures in two dimensional (2D) cell monolayers upon averaging over hundreds of division events. In addition, it is found that cell displacements exhibit super-diffusive behavior in a 2D confluent monolayer giavazzi2018tracking , where cell division and apoptosis occur continuously. By building on our previous studiesmalmi2018cell ; sinha2020spatially , which showed that imbalance between the rates of cell division and apoptosis, leads to unusual dynamics in a growing tumor spheroid, we investigate if a similar picture could be used to analyze the time-dependent cell trajectories generated in the 2D monolayergiavazzi2018tracking to generate new physical properties.

Here, we investigate the collective cell migration induced by division and apoptosis in a confluent 2D cell monolayer using experimental data on Madin-Darby Canine Kidney (MDCK) cellsgiavazzi2018tracking , in combination with an agent-based simulations malmi2018cell ; sinha2020spatially . We find that cell division leads to a directed motion of cells, which explains the super-diffusive behavior of cells with the mean-square displacement, $\Delta(t)\propto t^{1.5}$. The directed motion of the new born cells induce rapid rearrangements of the neighbors in the confluent tissue through topological (T1 and T2 transitions) changes. The cell rearrangement induced by a single cell division leads to vortex formation. With the occurrence of multiple cell births and deaths, a turbulent-like behavior naturally emerges with the mean size of the vortices being $\sim$ 100 $\mu m$. Interestingly, we find two scaling regimes for the kinetic energy spectrum, $E(k)$ with one following a Kolmogorov-like power law ($E(k)\approx k^{-5/3}$). The long-range correlated vortex motion leads to a long time tail (LTT) for the velocity auto-correlation function ($C_{v}(t)$), which decays as $C_{v}(t)\sim t^{-1/2}$, which differs from the well-known result $C_{v}(t)\sim\frac{1}{t\sqrt{lnt}}$ found in classical hard-disc-model fluids.

Vortex motion in a confluent cell monolayer: A confluent MDCK monolayer exhibits rich dynamics (see Figs. S1-S5 for some general features of the monolayer as it evolves in time) even though the cells are jammed. To probe the motility of the cells in the collective, we created the velocity vector field by applying the digital particle image velocimetrythielicke2014pivlab to the experimental datagiavazzi2018tracking . Collective motion of cells, in which a group of cells move in the same direction with similar magnitude of velocity (${\bf v}$)vicsek2012collective , occurs in different regions of the monolayer (see Fig. 1a). Interestingly, formation of several vortices are discernible in the velocity field, which indicates a coherent motion of cells (see the dashed rectangle in Fig. 1a as an example). We calculated the vorticity field (${\bm{\omega}}\equiv\nabla\times{\bf v}$), which describes the local rotational motion of the velocity flow (see Fig. 1b). Positive vorticity (yellow color) indicates clockwise rotation, while negative vorticity (blue color) shows anticlockwise rotation. The vortices in Fig. 1a are clearly elucidated in the vorticity field (see the bright yellow circle in the lower left of Fig. 1b and the dashed rectangle in Fig. 1a). Outward and inward velocity fluxes, characterized by divergence ($\text{div}({\bf v})\equiv\nabla\cdot{\bf v}$) of the velocity flow (Fig. 1c), are non-uniform in the monolayer.

Cell division as the driver of vortex motion: Previous theoretical and experimental studiesranft2010fluidization ; matoz2017cell ; malmi2018cell ; petridou2019fluidization have shown that active forces generated by cell division drives fluidization of tissues, leading to heterogeneous cell motility patterns both in a two-dimensional monolayer and in the three-dimensional spheroids. Such events (Fig. S1c) occur frequently in the experiments, as shown in Figs. 1a-c. Therefore, we surmised that the complex flow patterns observed here are driven by cell birth and death. To test this notion, we first show the positions of newly generated cells (born about 30 minutes before the snapshot in Fig. 1a) in the divergence field (see the purple dots in Fig. 1c). The origin of the outward flux (colored in yellow) is frequently accompanied by the appearance of newly created cells nearby. In addition, new cells arising from cell divisions are close to the vortices (see Figs. 1b-c). Therefore, it is likely that the flow patterns (Figs. 1a-c) in the monolayer are driven by active forces produced by the cell division.

Particle based simulations: To further elucidate the relevance of cell division, we simulated the dynamics of a cell monolayer using an agent-based model (details are in the Methods section and the SI). The velocity flow pattern in the simulations and experiments are similar (compare Figs. 1a and 1d). In addition, several pairs of vortices with opposite direction of rotation are also observed (see the lower middle of Fig. 1e). Note that in contrast to other studiesbelmonte2008self ; basan2013alignment ; bi2016motility , we do not include self-propulsion in our model. Thus, the observed velocity flow and vortex formation can only arise from self-generated active forces (SGAFs) arising from cell division. Indeed, cell division (see the purple dots in Fig. 1f), which produces an outward flow in the the divergence field (Fig. 1f), is close to the vortex location (see Figs. 1d-e). Hence, cell division and apoptosis generate complex cell motility patterns, resembling turbulent-like structures (Figs. 1a-f).

To illustrate how cell division leads to vortex formation, we tracked the collective cell motion in the dashed rectangle region in Fig. 1a for one hour (see Fig. 1g-l). The mother cell (in pink, see Fig. 1g at $t=100$ min) divides into two daughter cells (orange and blue) at $t=110$ min (see Fig. 1h). The vortex observed in Fig. 1a-b at $t=120$ min is established at $t=110$ min (see both the velocity and vorticity fields in Fig. 1h) immediately after cell division. The vortex persists for an additional 20 minutes (Fig. 1i-j) before vanishing (Fig. 1k-l). In contrast to previous experimental observationsrossen2014long ; doostmohammadi2015celebrating , where vortices are discernible only upon averaging over 100 division events, we found that a single cell division event leads to a vortex, lasting for more than 30 minutes. In addition, multiple cell division events drive a turbulent-like flow (Fig. 1a-f).

Directed motion of new born cells leads to anomalous diffusion: To study how vortices are formed by cell division explicitly, we started by monitoring the cell dynamics in the trajectories during the first $100$ minutes for the newly created cells (see Fig. 2a). For illustration purposes, we set the initial positions of all the cells at the origin. There is no obvious preferred direction of motion for cells, which suggests that they move isotropically (Fig. 2a). Surprisingly, the motion is directed after cell division over a distance of $\sim$ 10 $\mu m$ (see Fig. 2a and the inset figure for a few representative trajectories). The simulations also show a similar directed motion of cells during the first 100 minutes after cell divisions (see Fig. 2b).

To illustrate the consequence of cell-division induced directed motion, we calculated the mean-square displacement ($\Delta(t,t_{i})=\langle[\textbf{r}(t_{i}+t)-\textbf{r}(t_{i})]^{2}\rangle$) of cells at different initial times $t_{i}$, where r is the cell position. Both experiments and simulations (Figs. 2c-d) exhibit superdiffusive behavior, $\Delta(t,t_{i})\propto t^{\alpha}$, with an exponent $\alpha\approx 1.51\pm 0.1$ and $1.44\pm 0.1$ respectively, irrespective of the time windows used. The displacement distributions, ($P(\Delta x)$, $P(\Delta y)$), of cells calculated from both experiments and simulations (see Fig. S6a-d in the SI) deviate from the Gaussian distribution (see the green solid lines). It follows that cell division and apoptosis is the origin of superdiffusive motion of cells.

Cell division and apoptosis result in topological changes: The SGAFs sinha2020self due to cell division and apoptosis lead to a directed motion of newly generated cells, which influences the movement of cells nearby in a confluent tissue. Interestingly, there are several T1 transitions, which are critical during morphogenesisirvine1994cell , and cancer invasionlevayer2015cell , during a single cell division in the MDCK cell experiments. One example is shown in Figs. 2e-g (see also Fig. S7 in the SI). The force produced by the division of cell $1$ at $t=110$ min separates the two neighboring cells (cells 2 and 4) away from each other but brings cells $1$ and $3$ together at later times. Such rearrangement of cells $1-4$ and other cells nearby finally leads to vortex formation (see the velocity vector field in Figs. 2e-g and also Figs. 1g-l).

In addition to a single vortex induced by one cell division, (Figs. 1g-l), there is a pair of vortices, with opposite sense of rotation, in certain time frames (see the dashed rectangle in Fig. 3a). In the same region, there are three cell divisions, which lead to dramatic cell rearrangements nearby, 10 minutes before a vortex pair formation (the purple dots in the the dashed rectangle). To describe the formation of one vortex pair quantitatively, we plot the vorticity value along the arrow in the dashed rectangle in Fig. 3a using the digital particle image velocimetrythielicke2014pivlab . A smooth transition from anticlockwise to clockwise vorticity is clearly depicted (see Fig. 3b). Similarly, several pairs of vortices with opposite directions of rotation are found in simulations (see the lower middle of Fig. 1e). Besides, there is evidence for T2 transition (Fig. S8 in the SI), which arises from the extrusion (apoptosis) of a cell from the confluent monolayer. Taken together, these results suggest that cell division and apoptosis not only regulate cell numbers during morphogenesis and cancer development, they also drive dramatic cell rearrangements, as discovered recently for chick gastrulationfirmino2016cell , and germband extension of Drosophilada2007oriented . Finally, the cell rearrangement through topological T1/T2 transitions leads to the formation of vortices and turbulent-like structure in the confluent MDCK tissue.

Spatial correlations: To probe the spatial variations in the turbulent-like flow, we calculated space-dependent correlations associated using the velocity field, $C_{v}(R)$. From the cell position, $\textbf{r}_{i}(t)$, of the $i^{th}$ cell at different times, we calculated the velocity, $\textbf{v}_{i}(t)\equiv[\textbf{r}_{i}(t+\delta t/2)-\textbf{r}_{i}(t-\delta t/2)]/{\delta t}$, at time $t$ where $\delta t$ is the time interval. The spatial velocity correlation function ($C_{v}(R)$) is defined as,

$$C_{v}(R)=\langle\frac{\sum_{i}(\bf{v}(r_{i})-\langle v(r_{i})\rangle)\cdot(v(r_{i}+R)-\langle v(r_{i})\rangle)}{\sum_{i}(\bf{v}(r_{i})-\langle v(r_{i})\rangle)\cdot(v(r_{i})-\langle v(r_{i})\rangle)}\rangle.$$

(1)

The correlation function at different times (see Fig. 4a) has negative minimum at a distance $R_{v}\approx 100~{}\mu m$, which is about $10$ times greater than the mean cell size (see Fig. S3). The negative values for the correlation function are due to the anti-parallel velocity on the opposite sides of the vortices (see Fig. 1b). Therefore, $R_{v}$ is roughly the size of the vortices (Fig. 1), which are about 100 $\mu m$. The $C_{v}(R)$ calculated from our simulations supports such a link (see Fig. 4b). Interestingly, the collective motion of cells is quite robust, persisting for a long time, despite the increase in the cell density (Fig. S4a).

Long time tails (LTT): We calculated the velocity correlation function ($C^{i}_{v}(t)$) for cell $i$,

$$C^{i}_{v}(t)=\langle\frac{[\textbf{v}_{i}(t_{0}+t)-\langle\textbf{v}_{i}\rangle]\cdot[\textbf{v}_{i}(t_{0})-\langle\textbf{v}_{i}\rangle]}{[\textbf{v}_{i}(t_{0})-\langle\textbf{v}_{i}\rangle]\cdot[\textbf{v}_{i}(t_{0})-\langle\textbf{v}_{i}\rangle]}\rangle,$$

(2)

where $\langle\textbf{v}_{i}\rangle$ is the mean velocity of cell $i$, and the average is over time $t_{0}$ (24 hours). The ensemble average $C_{v}(t)$ for the monolayer is obtained by averaging $C^{i}_{v}(t)$ over all the cells. The decay of $C_{v}(t)$ follows $C_{v}(t)\sim t^{-\beta}$ with $\beta=0.4\pm 0.1$ (see the blue dash-dotted line in Fig. 4c). Our simulations also show a power-law decay with $\beta=0.59\pm 0.03$ (see Fig. 4d). Based on dimensional argument, we expect that the decay of $C_{v}(t)\sim\Delta(t)/t^{2}\sim t^{-0.5}$, which is fairly close to the experimental and simulation results. Using theory that accounts for SGAFs, the exponent $\beta=1/2$ in 2D (see the Methods section). A comparison from three different fit functions to the experimental and simulation data are also shown in the Fig. S9 in the SI. The emergence of LTT in $C_{v}(t)$ shows that motion of cells is persistent in a specific direction (Figs. 2a-d), which in turn is linked to vortex motion. To ascertain whether the power-law relation depends on the time interval selected to calculate $C_{v}(t)$, we varied $\delta t$. The value of $\beta$ stays around 1/2 (see Fig. S10 in the SI).

A long time tail relation is also found for the current-current correlation functions in disordered systems, which decays by $t^{-d/2}$ ( $d$ is the dimension of space)belitz1994anderson . Such a relation was first reported in simulations of classical hard-disc (sphere)-model fluidsalder1970decay , and was explained by correlated collision events theoreticallyernst1970asymptotic ; dorfman1970velocity . The difference in the value of $\beta$ between abiotic systems and MDCK monolayer is due to the generation of SGAFs due to cell division, which results in persistent (almost ballistic) directional motion of cells for long times.

Energy spectrum: It is clear that turbulent-like flow emerges naturally in the collective motion of cells (Figs. 1a-c). To characterize the nature of turbulent-like motion in the epithelial cells, we calculated the wave vector ($k$) dependent energy spectrum, $E(k)$, at different times (see Fig. 4e). There are two scaling regimes in $E(k)$ as a function of $k$. In the intermediate values of $k$, $E(k)\sim k^{-1.58(0.2)}$. The value of the exponent is close to the Kolmogorov-Kraichnan prediction, $-5/3$, found in the inertial turbulencekolmogorov1941local . We found a similar exponent value ($-1.4\pm 0.2$) from our simulations (see Fig. 4f). In the MDCK cells, the Reynolds number is small ($\approx 0$)marchetti2013hydrodynamics , which shows that the underlying mechanisms must be quite different for the turbulent-like motionalert2022active .

At smaller scales (large $k$), the results (Fig. 4e) obtained by analyzing the experimental trajectories show that $E(k)\sim k^{-\lambda}$ ($\lambda=3.5\pm 0.2$). From the agent-based simulations, we deduce that $E(k)\sim k^{-3.7(0.2)}$. Two comments are worth making: $(i)$ The $\lambda$ values from experiments and simulations are fairly close. It is interesting that scaling anticipated in the context of inertial turbulence is obeyed in living systems, given that energy injection, which is a consequence of cell division and apoptosis, is autonomous. $(ii)$ The prediction for inertial 2D turbulence at large $k$ is $\lambda\approx 3$, which differs from the estimated value for MDCK cells and simulations. It is unclear if the large $k$ finding in living tissues is simply a consequence of non-universal behavior or if large $k$ limit is not accessed in experiments and simulations. In active systems, such as microtubule-kinesin complex and bacterial suspensions, $\lambda$ ranges from 2 to 4.5wensink2012meso ; lin2021energetics ; guillamat2017taming ; alert2022active .

We investigated the dynamics of cell rearrangement and collective motion in a biologically relevant tissue friedl2009collective ; rorth2009collective ; poujade2007collective , by analyzing data from confluent MDCK cell monolayers, combined with agent-based simulations and theory. Cell division and apoptosis lead to a directed motion of cells for hours, which results in a universal super-diffusive behavior of cells characterized by the mean-square displacement exponent, $\alpha\approx 1.51$, irrespective of the initial time $t_{i}$ considered in the experimentsgiavazzi2018tracking . Our agent-based simulations yield $\alpha=1.44\pm 0.1$. It is worth emphasizing that, in the absence of cell division, there is complete cessation of motion in the simulations. Therefore, the complex dynamics in the simulations, which capture the salient features in the MDCK cells, arises solely due to self-generated active forces generated by cell division and apoptosismalmi2018cell ; sinha2020spatially ; sinha2020self .

The directed motion of newly generated cells induces rapid rearrangements of cells nearby through T1–cell intercalation and T2 transitions. Interestingly, cell rearrangements due to a single cell division result in vortex formation in confluent tissues, with lifetimes lasting for 30 minutes. More importantly, turbulent-like flux patterns emerges in the region where multiple divisions occur, which is also supported by our agent-based simulation model. From the spatial velocity correlation function, we estimate that the mean size of the vortices can be as large as to 100 $\mu m$, which is nearly 10 times larger than the size of a single cell ($\sim 10~{}\mu m$). Finally, the kinetic energy spectrum ($E(k)$) exhibits two distinct scaling behaviors of the collective cell movement. At intermediate values of $k$, $E(k)\sim k^{1.58\pm 0.2}$, which is close to the predicted Kolmogorov-Kraichnan behavior for the inertial turbulence. At large $k$ value, $E(k)\sim k^{-3.5}$, which does not seems to have counter part in inertial turbulence, and is likely to be non-universal. The exponent characterizing large $k$ behavior of $E(k)$ show substantial variations, depending on the systemalert2022active ; blanch2018turbulent ; mueller2019emergence .

In classical hard-disc-model fluidsalder1970decay ; ernst1970asymptotic ; dorfman1970velocity , a long time tail is expected for the velocity autocorrelation function with $C_{v}(t)\propto\frac{1}{t\sqrt{lnt}}$. Both experiments and simulations show a slower power law decay $C_{v}(t)\propto t^{-1/2}$ in the MDCK cell monolayer. To explain this finding, we developed a theory that accounts for active forces leading to the prediction of a slower decay, $C_{v}(t)\propto t^{-1/2}$. The theory and the agent-based simulations indicate that self-generated active forces lead to highly correlated responses. Given the simple model considered here, we expect that our studies can be used to describe collective motion of other cell types or a mixture of different cell types, and even other active systems in general.

In addition, we showed that the amplitude for the mean-square displacement (Fig. 2c) and the mean velocity of cells (Fig. S4b) decrease, as the initial time $t_{i}$ increases, which indicates a slowing down due to jamming or aging of the cell dynamics. Similar aging behavior has been noted in other types of cell monolayersgarcia2015physics . In addition, all the dynamical characteristics for cells, (superdiffusion, compressed exponential relaxation and aging) are found in many soft glassy systemsbouzid2017elastically ; gnan2019microscopic ; galloway2020scaling in the absence of active forces. Thus the physics used to describe abiotic glassy materials could be applied to understand living active matterkirkpatrick2015colloquium ; tjhung2020analogies . Finally, more detailed models by incorporating the size reductive divisionspuliafito2012collective , maturation and strengthening of cell-cell adhesiongarcia2015physics or other potential processes could be developed based on our present study to explain aging and other complex cell activities.

Summary of MDCK epithelial cell experimentsgiavazzi2018tracking : The MDCK (Madin-Darby Canine Kidney) epithelial cells were maintained in Dulbeccos Modified Eagle Medium, supplemented with $5\%$ fetal bovine serum and $1\%$ L-Glutamine. First, they were seeded in a 6-well plate to grow in complete medium until a confluent cell monolayer formed. Then, the images of the cell monolayer were acquired by Olympus IX81 inverted microscope (both in phase contrast and wide-field fluorescence) every minute for a 24-hour period. The phase-contrast microscopy was used to visualize the cell membrane and different organelles (see Fig. S1a). The nuclei of the cells at each time frame were also captured by wide-field fluorescence microscopy of EGFP-H2B expressing cells (see the cell nuclei shown in Fig. S1b for an example) and their positions were recorded from the single particle tracking by using the ImageJrueden2017imagej2 . A number of cell division events were observed in such a confluent monolayer (see Fig. S1c), while cell apoptosis/extrusion process led to the removal of cells from the monolayer. From the cell nuclei positions at each time frame, the cell trajectories can be obtained (see Fig. S1d, shows the first 150 minutes of recording for one of the field-of-views (FOV)). The cell trajectory in the figure exhibits a highly heterogeneous behavior. Some cells move in a rather straight and regular trajectory, while others exhibit irregular curly motion. The diverse movement is even clearer if the cell trajectories are plotted over the whole 24-hour period (see the two examples shown in the insets of Fig. S1d). The dynamically heterogeneous behavior is reminiscent of supercooled liquids, as was first shown by Angelini et al.angelini2011glass in the context of cell monolayers.

Particle-based simulations: We used a two-dimensional version of a particle-based model malmi2018cell ; sinha2020spatially to simulate cell dynamics. We mainly focused on investigating the roles of the cell division and apoptosis on the collective cell motion in order to provide microscopic insights into experimentsgiavazzi2018tracking . A single cell is modeled as a two-dimensional deformable disk. We use a periodic boundary conditions to mimic large systems. The cells interact with each other through repulsive (Hertzian contact mechanics) and attractive interactions (cell-cell adhesion). The motion of cells is described by an overdamped Langevin equation. The cells were allowed to grow and divide, as in the experiment. We also removed the cells from the monolayer at a constant rate to model the cell apoptosis and extrusion events that were observed in experiments. Because we do not consider self-propulsion, the cell motility can only arise from the self-generated active forces through cell growth, division and apoptosissinha2020self . The SI contains details of the model.

Theory of LTT driven by cell division: The turbulent-like motion in the MDCK cells in which active forces are generated by cell division, suggests that we use methods developed to treat turbulence in conventional fluids forster1977large ; dedominicis1979energy ; yakhot1986renormalization ; smith1992yakhot . The fundamental difference is that the random stirring that causes the turbulence in the conventional case is replaced by the internally generated cell division process in the biological case. Our aim is to compute the turbulent renormalized transport coefficients in the biological case, and from that determine the diffusive velocity autocorrelation function and the mean-squared-displacement. The relevant spatial dimension is $d=2$.

The model dynamical equations are,

$$\partial_{t}u_{i}+u_{j}\nabla_{j}u_{i}=-\nabla_{i}p+\nu_{0}\nabla^{2}u_{i}+f_{i}$$

(3)

$$\partial_{t}c+u_{j}\nabla_{j}c=D_{0}\nabla^{2}c+g$$

(4)

$$\nabla\cdot\mathbf{u}=0$$

(5)

Here, $u_{i}$ is the fluid velocity in the $i$-direction (repeated indices are to be summed over), $p$ is the pressure, $\nu_{0}$ is the bare kinematic viscosity, $c$ is the density of the active matter (MDCK cells), $D_{0}$ is a bare diffusion coefficient. The equations given above couple vorticity in the fluid motion to cell density. We have also imposed the incompressibility condition in Eq. (5). The statistical forces $f_{i}$ and $g$ are Gaussian noise terms that have zero mean and a nonzero two-point correlation functions that are local in space and time with gradient independent correlation strengths given by $\Delta_{f}$ and $\Delta_{g}$, respectively.

We treat the nonlinearities in these equations as a perturbation and work in Fourier space in both position and time, $(\mathbf{k},\omega)$. Further, the equations for both $u_{i}$ and $c$ are diffusive ¹¹1The pressure term in Eq. (3) is only used to enforce the incompressibility condition given by Eq. (5). and because $f_{i}$ and $g$ have identical statistical properties, we focus only on the mathematical structures and scaling solutions. We refer to the bare transport coefficients $D_{0}$, and their renormalized values by $D$. Structurally the one-loop contribution to $D$, denoted by $\delta D$, is given by,

$$\begin{split}&\delta D(\mathbf{k},\omega)\propto\Delta\int_{\mathbf{q},\omega_{1}}\frac{1}{[-i\omega-i\omega_{1}+D_{0}(\mathbf{k}+\mathbf{q})^{2}][\omega_{1}^{2}+D_{0}^{2}q^{4}]}\\ &\propto\Delta\int_{q>(\frac{\omega}{D_{0}})^{1/2},\omega_{1}}\frac{1}{[-i\omega_{1}+D_{0}q^{2}][\omega_{1}^{2}+D_{0}^{2}q^{4}]},\end{split}$$

(6)

where we have assumed that $\Delta\propto\Delta_{f}\propto\Delta_{g}$. In going from the first to the second line in Eq. (6) we have replaced the external frequency and wavenumber dependence in the integrand by a frequency cutoff. For $\mathbf{k}=0$, this can be justified in a scaling sense.

Carrying out the frequency and wave number integrals in $d=2$ gives,

$$\delta D(\omega\rightarrow 0)\propto\frac{\Delta}{D_{0}\omega}$$

(7)

Since $\delta D$ diverges as $\omega\rightarrow 0$, we use $\delta D\approx D$ and self-consistently replace $D_{0}$ in Eq. (7) by $D$. We conclude that in $d=2$,

$$D(\omega\rightarrow 0)\propto\frac{1}{\omega^{1/2}}.$$

(8)

This result in turn implies that the velocity autocorrelation function behaves as $C(t\rightarrow\infty)\propto 1/t^{1/2}$ and the mean-squared displacement is $\propto t^{3/2}$. These results are consistent with our numerical work as well as with experiments. In a similar manner, one obtains $D(\omega\rightarrow 0)\propto\frac{1}{\omega^{1/5}}$, which implies a velocity auto-correlation function $\propto 1/t^{4/5}$ in $d=3$.