Modeling Stochastic Variability in the Numbers of Surviving Salmonella enterica, Enterohemorrhagic Escherichia coli, and Listeria monocytogenes Cells at the Single-Cell Level in a Desiccated Environment

ABSTRACT Despite effective inactivation procedures, small numbers of bacterial cells may still remain in food samples. The risk that bacteria will survive these procedures has not been estimated precisely because deterministic models cannot be used to describe the uncertain behavior of bacterial populations. We used the Poisson distribution as a representative probability distribution to estimate the variability in bacterial numbers during the inactivation process. Strains of four serotypes of Salmonella enterica, three serotypes of enterohemorrhagic Escherichia coli, and one serotype of Listeria monocytogenes were evaluated for survival. We prepared bacterial cell numbers following a Poisson distribution (indicated by the parameter λ, which was equal to 2) and plated the cells in 96-well microplates, which were stored in a desiccated environment at 10% to 20% relative humidity and at 5, 15, and 25°C. The survival or death of the bacterial cells in each well was confirmed by adding tryptic soy broth as an enrichment culture. Changes in the Poisson distribution parameter during the inactivation process, which represent the variability in the numbers of surviving bacteria, were described by nonlinear regression with an exponential function based on a Weibull distribution. We also examined random changes in the number of surviving bacteria using a random number generator and computer simulations to determine whether the number of surviving bacteria followed a Poisson distribution during the bacterial death process by use of the Poisson process. For small initial cell numbers, more than 80% of the simulated distributions (λ = 2 or 10) followed a Poisson distribution. The results demonstrate that variability in the number of surviving bacteria can be described as a Poisson distribution by use of the model developed by use of the Poisson process. IMPORTANCE We developed a model to enable the quantitative assessment of bacterial survivors of inactivation procedures because the presence of even one bacterium can cause foodborne disease. The results demonstrate that the variability in the numbers of surviving bacteria was described as a Poisson distribution by use of the model developed by use of the Poisson process. Description of the number of surviving bacteria as a probability distribution rather than as the point estimates used in a deterministic approach can provide a more realistic estimation of risk. The probability model should be useful for estimating the quantitative risk of bacterial survival during inactivation.

T here is a risk of infection by foodborne pathogens when certain foods are ingested, even if as few as Ͻ10 to 100 Salmonella enterica or enterohemorrhagic Escherichia coli cells are present. Even one of these bacterial cells can cause foodborne disease, although the probability is very low (1). For example, a small number of cells in foods such as chocolate, salami, cheddar cheese, hamburger patties, and raw beef liver have been reported to present a medium risk of causing foodborne illness (2)(3)(4). Numerous predictive models based on deterministic approaches focusing on large bacterial populations, for example, more than 10 5 cells, have been developed to estimate the kinetics of inactivation of pathogenic bacteria (5). However, a deterministic approach results in limited predictions of bacterial behavior when dealing with low numbers of bacteria because this approach does not consider the variability and uncertainty of bacterial behavior (6). In a small population, the effect of the behavior of an individual bacterium becomes relatively large, and individual cell heterogeneity clearly appears when cell numbers are small (5). Because contamination of food with pathogens typically occurs with very low cell numbers, the use of probabilistic approaches that enable a description of the variability of the behavior of single cells is necessary to obtain more realistic estimates of the safety risk (7). Thus, there is a need to develop a predictive model to estimate the behavior of bacteria after the use of inactivation processes at the single-cell level.
In recent years, the need to consider the variability of the various factors that may influence predictive microbiology models has increasingly been recognized and has led to the development of more sophisticated stochastic models (8). Models that predict variability in the behavior of bacteria were developed by incorporating the probability distributions for variability or uncertainty model parameters in a Monte Carlo simulation (5,6,(9)(10)(11)(12). These models incorporated the variability caused by both the microorganism and the environment. Although bacterial behavior appears to fluctuate with low cell numbers, which may represent the natural stochastic variability in bacterial numbers, the randomness of the observed numbers of bacteria has not yet been directly considered or incorporated for evaluating bacterial behavior.
Recently, Koyama et al. (13) developed a sample preparation procedure for the probabilistic evaluation of bacterial behavior by obtaining bacterial numbers following a Poisson distribution (indicated by the parameter , which was equal to 2), which represents the variability in the occurrence of a natural event. In their paper, they suggested using the number of bacteria following a Poisson distribution ( ϭ 2) in a stochastic inactivation approach to investigate the variability in the numbers of surviving bacteria. In a similar approach, the Poisson distribution ( ϭ 2) was used to investigate the lag phase of single cells (14). To estimate the randomness of the number of bacterial cells that survive a process designed to kill bacteria, which could include heating, desiccation, or acid stress, we considered that the behavior of bacteria after the use of inactivation processes, incorporating the variability in the behavior of individual cells, can be described in a probabilistic model.
Foods with low water activity (a w ; i.e., foods with an a w of Ͻ0.85) do not support the growth of pathogenic bacteria. However, these foods could nevertheless be a vehicle for pathogenic bacteria and cause outbreaks of foodborne illnesses (2,15,16). Between 2007 and 2012, there were at least 7,315 cases of bacterial infection and 63 deaths worldwide due to the consumption of contaminated foods with a low a w (16), and the majority of these outbreaks were associated with S. enterica (2,15,17). Infection can result from the consumption of some low-a w foods, such as chocolate, salami, potato chips, and cheddar cheese, with as few as Ͻ10 to 100 S. enterica and enterohemorrhagic E. coli cells (2,4,(18)(19)(20). However, it is still difficult to predict bacterial survival at low levels with current experimental techniques. Estimation of the variability in the number of surviving bacteria will help to provide a more accurate estimate of the risk for contamination of low-a w foods.
The objective of the present study was to develop a probability model to estimate the variability in the number of surviving bacteria at the single-cell level in an environment with desiccation stress. In particular, we focused on describing the random-ness of the small numbers of bacteria that are present because of naturally occurring phenomena and that cannot be eradicated as variability in a Poisson process. We used bacterial cell numbers following a Poisson distribution ( ϭ 2) to determine the variability in the number of survivors at the single-cell level under desiccation stress at 5, 15, and 25°C. Furthermore, to confirm the validity of the experimentally observed naturally occurring randomness of the bacterial cell numbers in a Poisson process, we examined the bacterial death process using computer simulation. We determined whether a simulated number of survivors followed a Poisson distribution. The probability model developed in this study will enable estimation of the variability in the numbers of surviving bacteria and the probability of survival at the single-cell level.

RESULTS
Probability model for estimating the variability in the number of surviving bacteria over time in a desiccated environment. To describe the variability in the number of surviving bacteria with a Poisson distribution, we fitted a Poisson distribution to the bacterial number data during the bacterial death process using equation 2 (see Materials and Methods) and estimated the value of parameter of the Poisson distribution. Changes in the parameter of the Poisson distribution are described as an exponential function of the incubation time (Fig. 1). The root mean square error was below 0.095 under all conditions. There was no apparent dependency in the parameters, regardless of the bacterial strains examined (Table 1). Figure 2    a The parameter is calculated as follows: 10 Ϫbt n ϩ C , where b is the survival rate, t is the incubation time, n is the curvature of the kinetics of the Weibull distribution, and C is the initial bacterial distribution.
surviving bacteria at a specific time. The changes in the probability distribution of the numbers of surviving bacterial cells at the single-cell level during the bacterial death process were estimated as Poisson distributions ( Fig. 1 and 2). The model developed in the present study enabled us to estimate the changes in the distribution of bacterial numbers.
Computer simulation using random number generation for bacterial inactivation. We conducted a computer simulation to determine whether the frequency distribution of the 96 sets of bacterial numbers followed a Poisson distribution during the bacterial death process in a Poisson process (Fig. 3). Data for the 96 sets of bacterial cells following a Poisson distribution ( ϭ 2, 10, and 100) were examined as initial state {X N }, where N is the computer-generated total number of bacterial cells. X i represents the state space consisting of i surviving bacteria in a 96-well microplate (i ϭ N, N Ϫ 1, ..., 1). As shown in Fig. 3, a P value of Ͻ0.05 from the likelihood ratio test indicated that the hypothesis that the bacterial cell counts followed a Poisson distribution was rejected. We confirmed that among 1,000 simulations more than 82% of the simulated bacterial distributions ( ϭ 2) of each state {X N , X N Ϫ 1 , . . . , X 1 } followed a Poisson distribution during bacterial death. This result indicates that the distribution of surviving bacterial cells in a small population or at the single-cell level during the bacterial death process follows a Poisson distribution in the Poisson process. The simulated bacterial distribution ( ϭ 2) also followed a binomial or negative binomial distribution of more than 81%. These two distributions also describe the simulated bacterial distribution well. For the slightly high number of 10 cells, the simulation indicated that more than 80% of the cells followed a Poisson distribution, a binomial distribution, or a negative binomial distribution during the inactivation process. In contrast, for an initial number of 100 cells, less than 55% followed a Poisson distribution, a binomial distribution, or a negative binomial distribution. The mathematical description developed in this study will enable estimation of the probability of bacterial survival with small numbers of cells (e.g., 10 cells) along with its variability.

DISCUSSION
Previous studies have reported that in small populations, such as those comprising Ͻ100 cells, the behavior of individual cells can be very different from the overall behavior of the population (5,9). These previous studies have described variability from a biological perspective to be mainly the heterogeneity of individual cells, such as differences in tolerance to environmental stresses, such as heat, acidity, or desiccation. Here, we described variability to be the number of survivors based on mathematical descriptions. The focus of the present study was to describe the naturally occurring randomness of small numbers of bacteria on the basis of the assumption that the variability in the numbers of surviving bacteria followed a Poisson distribution. To assume a Poisson process, a Poisson distribution was used to describe variability in the number of surviving bacteria. As a result, the inevitable randomness of small numbers of surviving bacteria in the determination of the numbers was described successfully as the variability in the number of surviving bacteria.
We confirmed whether the distribution of survivors followed a Poisson, binomial, or negative binomial distribution. All three distributions fit well more than 80% of the simulated distributions with an initial state of 2 or 10 cells. Because the mean values for the targeted bacterial cell numbers were very small (2 or 10), a binomial distribution would be very similar to a Poisson distribution (21). In addition, when the sampling numbers for the negative binomial distribution were sufficiently increased (i.e., n ϭ 96), the relationship between the mean and variance of a negative binomial distribution is similar to that of a Poisson distribution (21). Thus, it can be concluded that the distribution of survivors is reasonably described as a Poisson distribution with small numbers of cells. The Poisson distribution is more useful than the other two distributions because the Poisson distribution is described by one parameter, .
The distinction in the behavior of a small number of bacterial cells compared with that of a larger bacterial population can be attributed to the law of large numbers; that is, as the sample size increases, the mean is more likely to approach the population average. However, to the best of our knowledge, the naturally occurring randomness of bacterial behavior that would represent variability in bacterial cell numbers has not been investigated in detail in a bacterial inactivation process to date. In this study, the variability was assumed to be the naturally occurring randomness in bacterial numbers. Given that the occurrence of such randomness cannot be determined as a fixed value,   The P values were estimated using a likelihood ratio test. If the P value was less than 0.05, the hypothesis that the bacterial cell counts followed a Poisson distribution was rejected.
it was described as a frequency distribution. Here, the variability in the number of surviving bacteria, which represents a naturally occurring random phenomenon, was described as a Poisson distribution. The model developed here enabled estimates of the randomness of bacterial numbers to reflect the realistic nature of random behavior. The growth behavior of bacterial cells at the single-cell level can be monitored via optical density measurements and/or under a microscope. To estimate the variability in cell growth, single-cell growth models have been developed (9,(22)(23)(24)(25)(26). However, there have been few studies in which these dynamics were monitored during the bacterial death process, owing to technical difficulties. Other than viable count methods, it is difficult to distinguish live cells from dead cells. One of the promising techniques for evaluating bacterial cell death processes is a time-lapse microscopic observation method (27). However, a special apparatus is needed for observation, and the procedure limits the applicable objective because the observation must take place while the cells are in a liquid environment. The sampling experiment examined in this study is technically easy to perform and more convenient to handle than monitoring of cell death. Therefore, to effectively and easily estimate changes in the probability of bacterial cell death during the bacterial death process, a sampling evaluation approach would be an effective alternative.
However, it is a challenge to obtain single cells, given the statistical variability and uncertainty associated with the experimental techniques (28,29). To overcome these technical difficulties with the collection of single cells, we modeled the bacterial numbers as a Poisson distribution, which is a natural distribution that can be easily reproduced (13). Thus, we described the number of surviving bacterial cells as a Poisson distribution rather than as a point estimate because the behavior of these survivors cannot be determined deterministically. Description of the numbers of surviving bacteria as a distribution rather than as the point estimates used in a deterministic approach can provide a more realistic estimate of risk.
Furthermore, we confirmed via computer simulation the variability in the number of surviving bacteria during the bacterial death process at the single-cell level (Fig. 3), which can be associated with the naturally random processes observed in many situations. The probability model developed here can be used to estimate the randomness of the bacterial numbers observed from the perspective of deterministically unexpected bacterial behavior. The variability in bacterial numbers, as defined in this study in the context of desiccation, and the concept of the random occurrence of surviving bacteria could be applied to other processes used to inactivate or kill bacteria, such as those involving heat, pressure, acidity, and irradiation. Probability models that describe the variability in the number of surviving bacteria as a probability distribution can be expected to play an important role in estimating the quantitative risk of bacterial survival. Listeria monocytogenes (ATCC 13932) was used. All the strains were maintained at Ϫ80°C in brain heart infusion broth (Merck, Darmstadt, Germany) containing 10% glycerol. A sterile metal loop was used to transfer the frozen bacterial cultures by scratching the surface of a frozen culture and placing it into 5 ml of tryptic soy broth (TSB; Merck, Darmstadt, Germany) in a sterile plastic tube. The cultures were incubated without agitation at 35°C for 24 h and were transferred, using a loop inoculum, at two successive 24-h intervals to obtain a more homogeneous and stable cell population. Grown cells were collected by centrifugation (1,000 ϫ g, 10 min at 25°C), and the resulting pellet was washed twice with peptone water and subsequently resuspended in 5 ml of peptone water.

MATERIALS AND METHODS
Preparation of single cells. The procedure used for the preparation single bacterial cells was that of Koyama et al. (13). Briefly, single cells of each pathogenic bacterium were prepared by the use of 10-fold dilutions of bacterial cultures. The initial cell count was assumed to be approximately 10 9 CFU/ml. The inoculum was further diluted (10-fold dilution series) in peptone water to obtain a solution with 10 3 CFU/ml, and the bacterial cell number was determined via direct plating of 100 l of the inoculum of a sample with 10 3 CFU/ml onto five tryptic soy agar plates, which were then incubated at 35°C for 24 h. The data are expressed as the mean value (x ) for five replicates of plate counts. An aliquot of the solution with 10 4 CFU/ml stored from the previous day of the experiment was diluted x /10-fold in pure water to obtain a solution with 1,000 CFU/ml. Pure water was used to prevent the bacteria from being stressed by the high osmotic pressure used during desiccation. This procedure resulted in a final solution of 2 CFU/2 l.
Assumption of a stochastic bacterial death process. In the present study, we assumed that the bacterial death process followed a Poisson process. An occurrence or event was independent of itself, and occurrences or events appeared one by one over an infinite time period. The variability in the total number of occurrences or events until a specific time period follows a Poisson distribution in a Poisson process (21,30). According to the nature of a Poisson process, we experimentally investigated and computationally simulated bacterial death.
Bacterial inactivation in a desiccated environment. Single cells of each bacterial serotype in 2 l of solution were dispensed into each well of 96-well microplates. The prepared microplates were dried in a drying chamber using silica gel (10% to 20% relative humidity) at 5, 15, and 25°C. We confirmed the survival or death of the cells in each well by adding 100 l of TSB about every 1 or 2 h until 10 h and every 24 h until 48 h. One microplate was used to determine the survival or death of cells at each time and under each condition. The microplates to which TSB was added were stored at 25°C for 1 week to allow the recovery of injured bacterial cells (31,32). The turbidity of the wells was observed to determine the survival or death of cells.
Analysis of data on the change in the survival rate during the drying process. We assumed that the number of survivors followed a Poisson distribution during bacterial death. We determined the proportion of wells that contained surviving cells by counting the number of wells in the 96-well plates with surviving bacteria to estimate the survival probability parameter of the Poisson distribution. The Poisson distribution is described by equation 1: where P(X ϭ k) is the probability of finding k cells in a randomly chosen well and is the parameter of the Poisson distribution. From the rate of survival determined from wells with surviving bacteria, we estimated the parameter of the Poisson distribution. From equation 1, the parameter of the Poisson distribution was calculated as follows: Given that the parameter of the Poisson distribution was continuous and nonnegative, we used a nonnegative decreasing mathematical function to analyze the data. We described a change in the value of parameter of the Poisson distribution over time as an exponential function on the basis of the Weibull distribution by nonlinear regression. The Weibull distribution is commonly used to describe the log reduction (33) and has been used to describe bacterial inactivation kinetics (34,35). The mean value of the response variable , which is the number of survivors, follows a Poisson distribution of mean . Equation 3 was used as the basic fitting model: where b and n represent the survival rate and the curvature of kinetics of the Weibull distribution, respectively, and C and t are the initial bacterial distribution and the incubation time, respectively, and were estimated by use of the parameter. We used equation 3 to estimate the parameter of the Poisson distribution as follows: The variability in the observed outcomes was estimated as the variability in the mean value for the 96 sets of bacterial numbers generated from the Poisson distribution. We estimated the variability in the observed outcomes as the 95% confidence interval of the mean value generated from the Poisson distribution.
Bacterial inactivation simulation using random number generation. We simulated the stochastic death process for 96 sets of bacterial cells by the generation of random numbers by computer simulation in a Poisson process (Fig. 4). We determined whether the frequency distribution of 96 sets of bacterial cell numbers followed a Poisson distribution during the bacterial death process. The simulation of the bacterial death process was based on the following presuppositions: the time to inactivation of an individual cell is independent of any surrounding cells, and the time to inactivation of individual cells occurs one by one over an infinite time period.
The computer simulation examined a bacterial death process from 96 sets of virtual bacterial cells following a Poisson distribution. Ninety-six random numbers following a Poisson distribution ( ϭ 2) were generated, as shown in Fig. 4a and 5a. {X N } represents the initial state space consisting of the number of bacteria in a 96-well microplate. N is the computer-generated total number of bacterial cells. We uniformly and randomly chose a surviving bacterial cell from {X N } and assumed it to be a dead cell. The number of survivors then became N Ϫ 1, and the initial state {X N } transitioned to state {X N Ϫ 1 } (Fig. 4b). As the number of survivors decreased one by one, transition of the state occurred until state {X i } (where i ϭ N, N Ϫ 1, . . ., 1) ( Fig. 4c and 5b) became state {X 1 }, which contained one surviving cell ( Fig. 4d and 5c). The computer simulation was repeated 1,000 times. The algorithm is summarized as a flowchart in Fig. 5. We repeated the computer simulation and examined the bacterial death process from 96 sets of virtual bacterial cells following two different Poisson distributions ( ϭ 10 or ϭ 100).
We determined whether the probability distribution of surviving bacterial cell numbers in each state {X N , X N Ϫ 1 , . . . , X 1 } followed a Poisson distribution during the death process using a likelihood ratio test (36). If the P value was less than 0.05, the hypothesis that the bacterial cell counts followed a Poisson distribution was rejected. We also determined whether a binomial distribution and a negative binomial distribution fit the probability distribution of surviving bacterial cell numbers in each state.