Electrical signaling in three-dimensional bacterial biofilms using an agent-based fire-diffuse-fire model

Agent-based models were used to describe electrical signaling in bacterial biofilms in three dimensions. Specifically, wavefronts of potassium ions in Escherichia coli biofilms subjected to stress from blue light were modeled from experimental data. Electrical signaling occurs only when the biofilms grow beyond a threshold size, which we have shown to vary with the ${\mathrm{K}}^{+}$ ion diffusivity, and the ${\mathrm{K}}^{+}$ ion threshold concentration, which triggered firing in the fire-diffuse-fire model. The transport of the propagating wavefronts shows superdiffusive scaling on time. ${\mathrm{K}}^{+}$ ion diffusivity is the main factor that affects the wavefront velocity. The ${\mathrm{K}}^{+}$ ion diffusivity and the firing threshold also affect the anomalous exponent for the propagation of the wavefront determining whether the wavefront is subdiffusive or superdiffusive. The geometry of the biofilm and its relation to the mean-square displacement (MSD) of the wavefront as a function of time was investigated for spherical, cylindrical, cubical, and mushroom-like structures. The MSD varied significantly with geometry; an additional regime to the kinetics occurred when the potassium wavefront leaves the biofilm. Adding cylindrical defects to the biofilm, which are known to occur in E. coli biofilms, also gave an extra kinetic regime to the wavefront MSD for the propagation through the defect.

Figure 1

Schematic of the BSim agent-based model. The two main components of the model are shown: the agents and the process by which agents interact. The agents are shown as green spheres, and the chemical field is shown in blue. BSim allows the partial differential equations for the diffusion of the chemical field to be calculated [21].

Figure 2

An example of a spherical biofilm constructed using BSim. The bacteria are stationary and are shown in blue. The EPS is neglected in the model. The bacteria are randomly placed within the sphere, but bacterial overlaps are not permitted. The simulation box is a cube (gray).

Figure 3

(a) Simulation of the release of potassium ions of a single bacterium over different time intervals. A single spike occurs in the middle, where the single cell was located as expected from Fick's laws of diffusion, and a decrease in the peak heights occurs over time. A Gaussian curve of the form $\alpha exp[-{((x-\beta )/\gamma )}^2]$ was fit for each of the time curves. For 100, 200, 300, 400, and 500 s the values for $\alpha$ were $3.33\times {10}^6, 1.10\times {10}^6, 5.66\times {10}^5, 3.50\times {10}^5$, and $2.39\times {10}^5$, respectively; the values for $\beta$ were 2.54, 2.54, 2.55, 2.55, and 2.56, respectively; the values for $\gamma$ were 9.64, 13.6, 16.6, 19.1, and 21.3, respectively. The fits all corresponded to a $R^2$ value of 0.9999. (b) The logarithm of the intensity and the square of the distance were considered to ensure it is Gaussian. The straight lines suggest that the trend is of Gaussian nature, though it is slightly off-center due to the discrete nature of the lattice, i.e., a one-pixel bias. The potassium diffusivity ($D_c$) was $0.1\mu {\text{m}}^2/\text{s}$ and the decay constant ($k_s$) was set to a value of $4.44\times {10}^{-4}$. The quantity of potassium ions added ($\sigma$) was $5\times {10}^9$.

Figure 4

Average velocity of a wavefront of ${\mathrm{K}}^{+}$ ions as a function of ion diffusivity, averaged over three simulations. The decay rate ($k_s$) was kept constant at $4.44\times {10}^{-4}$ molecules per second and the ${\mathrm{K}}^{+}$ ions released was kept constant at $5\times {10}^9$ molecules in the FDF model for the whole biofilm. The firing threshold ($c^{*}$) was maintained at $7\times {10}^6$. The method for determining the intensities was the spherical shell method. The error bars show the standard error on the mean of three runs. The power-law fit was of the form $f\left(x\right)=a{x}^b+c$, where $a=22.3, b=1.85$, and $c=0.75$. It had a regression ${\mathrm{\Lambda }}^2$ value of 0.9831.

Figure 5

The instantaneous velocity of ${\mathrm{K}}^{+}$ ion wavefronts as a function of ${\mathrm{K}}^{+}$ diffusivity calculated for a range of radii. An average of three simulations is shown. The decay rate ($k_s$) was kept constant at $4.44\times {10}^{-4}$ molecules per second, and the amount of potassium ions released ($\sigma$) was kept constant at $5\times {10}^9$ molecules in the FDF model. The firing threshold ($c^{*}$) was maintained at $7\times {10}^6$. The error bars are the standard error on the mean of three simulations. The power-law fit of the average for each diffusivity point was of the form $v\left(D\right)=a{D}^b+c$, where $a=4.21, b=0.79$, and $c=0.014$. It had a regression ${\mathrm{\Lambda }}^2$ value of 0.9958.

Figure 6

Power-law scaling exponents of the instantaneous velocity of the ${\mathrm{K}}^{+}$ ion wavefront as a function of radial position inside the biofilm for an average of three runs (the fractional radius = 1 on the edge of the biofilm). The decay rate ($k_s$) was kept constant at $4.44\times {10}^{-4}$ molecules per second, and the ions released ($\sigma$) were kept constant at $5\times {10}^9$ in the FDF model. The error bars are the standard error mean of the three runs. The power law was in the form of $f\left(R\right)=a{R}^b+c$, where $a=0.97, b=0.65$, and $c=0.23$. It had a regression ${\mathrm{\Lambda }}^2$ value of 0.9987.

Figure 7

Globally averaged fluorescence intensity from (a) experimental data compared to (b) simulation data. (a) ThT was used as a Nernstian probe in fluorescence microscopy experiments with an E. coli biofilm [6]. The error bars in the top subfigure represent the standard deviation over more than three experiments. (b) The global oscillations for the bottom figure are taken for a simulation with ${\mathrm{K}}^{+}$ diffusivity ($D_c$) of $0.1\mu {\mathrm{m}}^2$, decay rate ($k_c$) of $5\times {10}^{-3}$ molecules/s, $5\times {10}^9$ potassium ions ($\sigma$) added per spike, and $1\times {10}^7$ potassium ions as the firing threshold ($c^{*}$) for the FDF mechanism. The first wavefront of the FDF stops at $t_1=180$ sec and the second wavefront begins (a new spike) at $t_2=1440$ sec. The simulation stops at $t_3=3600$ sec. 

Figure 8

Mean-square displacement of (a) centrifugal and (b) centripetal ${\mathrm{K}}^{+}$ wavefronts for a simulation with ${\mathrm{K}}^{+}$ diffusivity of $0.1\mu {\mathrm{m}}^2$, decay rate ($k_s$) of $5\times {10}^{-3}$ molecules/s, $5\times {10}^9$ quantity of ${\mathrm{K}}^{+}$ ($\sigma$) added per hyperpolarization event and $3\times {10}^9$ molecules as the threshold ($c^{*}$) for the firing mechanism. (a) The curve fit is Eq. (4), where $\gamma$ is $1.21\pm 0.12$ and $R_c$ is $6.17\pm 1.84\mu \text{m}$. The regression $R^2$ value is 0.9911. The residuals of the model are plotted underneath the main graph. (b) The curve fit is again Eq. (4), where $\gamma$ is $2.26\pm 0.31$ and $R_c$ is undefined due to the nature of the centripetal wavefront. The regression ${\mathrm{\Lambda }}^2$ value is 0.9840. The residuals of the model are plotted underneath the main graph. (c) Logarithmic plot for the time and $R^2\text{−}{c}^2$ for the centrifugal wavefront; (d) logarithmic plot for the time and $R^2\text{−}{c}^2$ for the centripetal wavefront.

Figure 9

Instantaneous velocity as a function of time of (a) the centrifugal and (b) centripetal ${\mathrm{K}}^{+}$ wavefronts for simulations where the ${\mathrm{K}}^{+}$ diffusivity ($D_c$) is $0.1\mu {\mathrm{m}}^2$, the decay rate ($k_c$) is $5\times {10}^{-3}$ molecules/s, the quantity of ${\mathrm{K}}^{+}$ added per hyperpolarization event ($\sigma$) is $5\times {10}^9$, and $3\times {10}^9$ molecules as the firing threshold ($c^{*}$) for the FDF mechanism.

Figure 10

Instantaneous velocity as a function of curvature ($=1/R$) for (a) the centrifugal and (b) the centripetal wavefronts from simulations with ${\mathrm{K}}^{+}$ diffusivity ($D_c$) of $0.1\mu {\mathrm{m}}^2$, decay rate ($k_s$) of $5\times {10}^{-3}$ molecules/s, $5\times {10}^9$ quantity of ${\mathrm{K}}^{+}$ ($\sigma$) added per hyperpolarization event, and $3\times {10}^9$ molecules as the firing threshold ($c^{*}$) for the FDF mechanism. The radius is defined as zero at the center of the biofilm, and the biofilm periphery is the maximum radius. The error bars show the 95% confidence interval for the velocities.

Figure 11

Mean-square displacement of (a) the centrifugal and (b) centripetal ${\mathrm{K}}^{+}$ wavefronts as a function of time measured using confocal microscopy with ThT fluorophores. (a) The curve fit is Eq. (4), where $\gamma$ is $1.22\pm 0.15$ and $R_c$ is $4.71\pm 0.98\mu \mathrm{m}$. (b) The curve fit is Eq. (4), where $\gamma$ is $2.43\pm 0.08$ and $R_c$ is undefined due to the nature of the centripetal wavefront.

Figure 12

Instantaneous velocity as a function of radius of (a) the centrifugal and (b) the centripetal ${\mathrm{K}}^{+}$ wavefronts measured from the first peak of the confocal microscopy data. An inverse relation of the ${\mathrm{K}}^{+}$ wavefront velocity is observed with radius for (a), whereas (b) has a minimum. The error bars are the standard deviation of more than three experiments.

Figure 13

Instantaneous velocity as a function of time for (a) the centrifugal and (b) the centripetal ${\mathrm{K}}^{+}$ wavefronts measured from the first hyperpolarization event of the confocal microscopy data. An inverse relationship of the wavefront velocity with time is shown for (a), whereas (b) has a minimum. The error bars are the standard deviation of more than three experiments.

Figure 14

(a) Mean-square displacement of ${\mathrm{K}}^{+}$ wavefronts as a function of time for different ${\mathrm{K}}^{+}$ diffusivities from ABM simulations. The different colors represent the different ${\mathrm{K}}^{+}$ diffusivities shown in the legend. All ${\mathrm{K}}^{+}$ diffusivities are given in $\mu m^2/s$. The dashed trendlines for different $\gamma$ power-law scaling exponents are also shown for reference. (b) The plot in (a) but scaled in a logarithmic manner (base 10).

Figure 15

Mean-square displacement of the ${\mathrm{K}}^{+}$ wavefronts as a function of time for different ${\mathrm{K}}^{+}$ concentration firing thresholds ($c^{*}$) from ABM simulations. The different colors represent the different firing thresholds shown in the legend. All firing thresholds are given in terms of the number of molecules needed to trigger the FDF mechanism.

Figure 16

Critical radius for propagation of ${\mathrm{K}}^{+}$ wavefronts as a function of the diffusion coefficients of the ${\mathrm{K}}^{+}$ ions ($D_c$, the diffusivity) from ABM simulations. An average of three runs is shown. The critical radii decrease as the ${\mathrm{K}}^{+}$ diffusivity increases. The error bars are the standard error mean of the three runs.

Figure 17

Different bacterial biofilm geometries that were used for the FDF ABM in addition to a spherical geometry (Fig. 2): (a) a cylinder, (b) a cube, and (c) a mushroom shape.

Figure 18

Total intensity of the potassium ions averaged over an entire biofilm as a function of time for different biofilm shapes spiked at their centers. The parameters were kept constant throughout the different geometries: the ${\mathrm{K}}^{+}$ ion diffusivity ($D_c$) was $0.1\mu {\mathrm{m}}^2/s$, the total potassium ions added ($\sigma$) was $5\times {10}^9$, the firing threshold ($c^{*}$) was $1\times {10}^7$, and the decay rate ($k_s$) was $5\times {10}^{-3}$. The first wavefront of the FDF stopped at $t_1=180$ s, and the second wavefront began (new hyperpolarization event) at $t_2=1440$ s. The simulation stopped at $t_3=3600$ s.

Figure 19

Schematic diagram of a cylindrical biofilm with a cylindrical defect shown in red used in the simulation. The bacteria are randomly placed in the cylinder except for the area shown in red. (a) Top view of the plane perpendicular to the $y$ axis. (b) Side view of the plane parallel to the $y$ axis. The blue arrows indicate the directions in which the wavefront propagation was quantified.

Figure 20

Mean-square displacement, $R^2$, of the ${\mathrm{K}}^{+}$ wavefront as a function of time for (a) a cylinder, (b) a cube, (c) a mushroom-like shape, and (d) a cylinder with a cylindrical defect. The simulation parameters were kept constant for the different geometries. The ${\mathrm{K}}^{+}$ diffusivity ($D_c$) was $0.1\mu {\mathrm{m}}^2/s$, the total ${\mathrm{K}}^{+}$ added ($\sigma$) was $5\times {10}^9$, the firing threshold ($c^{*}$) was $1\times {10}^7$, and the decay rate ($k_s$) was $5\times {10}^{-3}$. The trend lines for different $\gamma$ values are shown for the different regimes. The error bars are the standard mean error on the wavefronts for three runs.

Figure 21

Mean-square displacement, $R^2$, of the ${\mathrm{K}}^{+}$ wavefront as a function of time for (a) a small cylinder, (b) a small cube, and (c) a small cylinder with a cylindrical defect. The simulation parameters were kept constant for the different geometries. The ${\mathrm{K}}^{+}$ diffusivity ($D_c$) was $0.1\mu {\mathrm{m}}^2/s$, the total ${\mathrm{K}}^{+}$ added ($\sigma$) was $5\times {10}^9$, the firing threshold ($c^{*}$) was $1\times {10}^7$, and the decay rate ($k_s$) was $5\times {10}^{-3}$. The error bars are the standard mean error on the wavefronts for three runs.

Figure 22

Concentration profile of potassium ions released from a single bacterium as a function of time ($t$). ${\rho }_c$ is the potassium concentration density, i.e., the ${\mathrm{K}}^{+}$ concentration per unit time. The area of the rectangles equals the total ${\mathrm{K}}^{+}$ concentration released. (a) The initial FDF model with instantaneous release and a finite time step. (b) A sustained release time ($r$) and a dormancy time ($d$) are shown.

Figure 23

Mean-square position ($R^2$) of the ${\mathrm{K}}^{+}$ wavefront, $R^2$, as a function of time for different rise and dormancy times from ABM simulations. The ${\mathrm{K}}^{+}$ diffusivity ($D_c$) was $0.1\mu {\mathrm{m}}^2/s$, the total ${\mathrm{K}}^{+}$ added ($\sigma$) was $5\times {10}^9$, the firing threshold ($c^{*}$) was $1\times {10}^7$, and the decay rate ($k_s$) was $5\times {10}^{-3}$. The trend line is a power law with exponent $\gamma$ = 2 (dashed line) showing ballistic scaling.