Abstract
For manufacturers of sterile drug products, steam sterilization is a common method used to provide assurance of the sterility of manufacturing equipment and products. The validation of sterilization processes is a regulatory requirement and relies upon the estimation of key resistance parameters of microorganisms. Traditional methods have relied upon point estimates for the resistance parameters. In this paper, we propose a Bayesian method for estimation of the well-known DT, z, and Fo values that are used in the development and validation of sterilization processes. A Bayesian approach allows the uncertainty about these values to be modeled using probability distributions, thereby providing a fully risk-based approach to measures of sterility assurance. An example is given using the survivor curve and fraction negative methods for estimation of resistance parameters, and we present a means by which a probabilistic conclusion can be made regarding the ability of a process to achieve a specified sterility criterion.
LAY ABSTRACT: For manufacturers of sterile drug products, steam sterilization is a common method used to provide assurance of the sterility of manufacturing equipment and products. The validation of sterilization processes is a regulatory requirement and relies upon the estimation of key resistance parameters of microorganisms. Traditional methods have relied upon point estimates for the resistance parameters. In this paper, we propose a Bayesian method for estimation of the critical process parameters that are evaluated in the development and validation of sterilization processes. A Bayesian approach allows the uncertainty about these parameters to be modeled using probability distributions, thereby providing a fully risk-based approach to measures of sterility assurance. An example is given using the survivor curve and fraction negative methods for estimation of resistance parameters, and we present a means by which a probabilistic conclusion can be made regarding the ability of a process to achieve a specified sterility criterion.
1. Introduction
In the pharmaceutical industry, moist heat (or steam) is a commonly used sterilizing agent for thermal destruction of microorganisms. Steam sterilization processes are used to claim sterility of manufacturing equipment such as aseptic filler piping used for a sterile drug product. The sterilization process is governed by time, temperature, and pressure, and it is typically implemented in three phases: heating; sterilizing, and cooling. Validation of a sterilization process can rely upon a biological and/or mathematical (physical) approach (1). The biological and mathematical approaches for measuring the efficacy of the sterilization process are complementary to each other (2). Furthermore, it is typically recommended that lethality be determined using both methods so that they “validate” each other. In this paper, we propose a Bayesian approach to the physical method of assessing the ability of a sterilization process to achieve a specified sterility criterion.
Sterility assurance is fundamentally a probability measure. It is defined as the probability of a single viable microorganism occurring on an item after sterilization, or the probability of a non-sterile unit (PNSU) (1). In the pharmaceutical industry, a system is considered to be “sterile” when PNSU ≤10−6. Note that the definition of sterile does not imply that one in a million products is allowed to be non-sterile, but it admits a finite statistical probability that a microorganism may survive the sterilizing process (2, 3). Because probability is a representation of uncertainty, the application of Bayesian methods to the field of sterilization validation makes sense. A Bayesian approach requires that all uncertainty be represented by probability. Traditional sterility measures, such as the well-known decimal reduction (DT) value and the viable count of a biological indicator (BI), are estimated with a high level of uncertainty due to many factors relating to process parameters that are subject to kinetic and biological variations. In discussing methods for gathering and analyzing microbial heat-destruction data, Pflug et al. (4) note that, in the biologic world, there are no “standard” microorganisms and that biologic systems are well known for their variability. Recovery or growth mediums used for enumerating microorganisms are also biologic products with variability from lot to lot, and sampling techniques tend to be susceptible to environmental variations. Although the uncertainties in microbiological parameters have been noted for some time in the literature, inclusion of variability in thermal process calculations has been reported infrequently (5). For these reasons, large margins of error are permitted for these experimentally determined parameters. The ISO 11138-1 standard (6), for example, allows for the tested DT value of a BI to be within ± 20% of the manufacturer's stated value and the viable count of the BI to be between 50% and 300% of the stated value. However, this uncertainty is not taken into account in final sterility assurance calculations.
In this paper, we present a Bayesian approach to the estimation of DT, z, and Fo values that are commonly used in the pharmaceutical industry to provide evidence of sterility assurance. Each of these parameters will be defined and their theoretical derivations given in the following section. We show in this paper that by using Bayesian methods to quantify uncertainty in sterility calculations, a truly probabilistic or “risk-based” conclusion regarding the effectiveness of a sterilization process can be made. Moreover, we show that a Bayesian approach provides a coherent framework for combining knowledge gained about DT from multiple estimation procedures into a single probability distribution used for inference. We begin in Section 2 by providing the theoretical framework for steam sterilization processes. Section 3 presents a brief introduction to Bayesian estimation concepts, and in Section 4 we show how Bayesian methods can be used to provide probability-based estimates and intervals for the parameter values DT, z, and Fo. In Section 5 we define the risk of non-sterility in terms of the PNSU measure and present a framework for quantifying such risk a posteriori. A complete example is given in Section 6, where the survival curve method and the Holcomb Spearman-Karber procedure are used to estimate DT and a sterilization process is evaluated using the posterior probability distribution of Fo. Finally, some concluding remarks are given in Section 7.
2. Theoretical Background
The inactivation of microorganisms exposed to constant thermal stress follows the differential equation for a first-order reaction or birth-death process (see 1 and 7) given by
where N(t) is the number of microorganisms at time t, and k is the inactivation constant for the microorganism. Although k is typically not directly estimated, it is dependent on temperature in steam sterilization via the Arrhenius law (see 8, 9),
where δ denotes the constant “pre-exponential” factor, E is the activation energy (J mol–1), τ is the temperature (K), and R denotes the universal gas constant [8.31 J/(mol K)]. From eq 1, it follows that
where N(0) denotes the initial number of organisms and k′ = k/ln(10). The decimal reduction time, DT, is defined as the time, at a given temperature, T (°C), required to reduce the number of organisms by a factor of 10. Because N(DT) = N(0)/10, it follows from eq 2 that DT = 1/k′ and
where log denotes logarithm of base 10. Equation 3 shows that the logarithm of the number of surviving organisms is linearly related to time. This is the theoretical survival curve for the organism.
A sterility criterion, S, is typically defined in terms of log reduction. For example, a sterilization process that achieves a 6 log reduction reduces the initial population of organisms by a factor of 1,000,000. As a result, S = log[N(0)/N(t)] and the minimum time, FT, to satisfy the sterility criterion at a given temperature, T, can be calculated using eq 3 as
The resistance coefficient, z, is the temperature (°C) increase required to reduce DT by a factor of 10. It is assumed that T and logDT are linearly related and that the z value is the negative reciprocal of the slope of the line, which is often referred to as the thermal resistance curve (1). Thus, if the reference decimal reduction time is DR at temperature TR, then
or
As steam sterilization processes are typically carried out within a small temperature range, the z value is usually considered constant for practical purposes (1). Following from eq 4, FT at resistance coefficient z is
Lethality is a measure of killing effect. For comparison purposes, a unit of lethality is defined as the killing effect of 1 min (FR = 1) at a reference temperature, TR, usually taken to be 121.1 °C. Thus, if T = 111 °C and z = 12 °C, then the killing effect, F11112, is approximately 0.1 of equivalent minutes at TR = 121.1 °C. The term 10(T−121.1)/z in eq 6 is typically defined as the lethal rate.
As sterilization processes have heat-up and cool-down phases, they do not reach the reference sterilization temperature (121.1 °C) instantaneously. As result, the inactivation of microorganisms is a continuous process during which the rate of inactivation (DT value) changes with temperature. It is therefore possible to “accumulate” lethality over the duration of the sterilization process. The accumulated lethality, with respect to the reference temperature, is denoted in the pharmaceutical industry by Fo and is calculated as
where n is the number of samples, Ts is the temperature recorded at sample s, and Δt denotes the sampling interval in minutes.
Before presenting Bayesian methods for estimation of DT, z, and Fo, we first give a brief overview of the Bayesian approach to statistical inference in the following section.
3. Bayesian Estimation
Bayesian estimation uses probability to quantify uncertainty. Suppose θ is a parameter of interest. In Bayesian methods, a prior distribution, π(θ), represents all that is known prior to the current experiment. The probability distribution of the data, y, is written as a function of θ called the likelihood, denoted by l(θ|y). Bayes' theorem combines these two sources of information, updating the prior to form the posterior probability distribution
where the integral is over the possible values of θ, here assumed to be continuous. Of course, for models with more than one parameter, the integral in the denominator becomes a multiple integral. Due to the intractability of the integrals in most “real-life” problems, integration is almost always done using iterative simulation methods. Specialty software such as OpenBUGS (www.openbugs.net), JAGS (http://mcmc-jags.sourceforge.net/), and Stan (http://mc-stan.org/) can be used in conjunction with R (https://www.r-project.org/) to perform Bayesian inference in many types of models. The SAS procedure PROC MCMC is also useful. The posterior distribution represents the state of knowledge after the experiment is performed and can be used as a prior for a subsequent experiment. The data therefore “update” the current state of knowledge after each experiment. As Gelman et al. (10) note, a primary motivation for Bayesian thinking is that it facilitates a common sense interpretation of statistical conclusions. For example, a Bayesian interval estimate for an unknown parameter of interest (known as a credible interval) can be directly regarded as having a high probability of containing the unknown parameter value. Quantities such as posterior moments, credible intervals, and tail probabilities are generally used to perform Bayesian statistical inference.
Examples of the application of Bayesian methods in pharmaceutical science and technology can be found in several PDA Journal articles—see, for example, Yang and Zhang (11), Yang (12), and Claycamp (13). In exploring the use of probability concepts in pharmaceutical quality risk management, the U.S. Food and Drug Administration's Claycamp (13) suggests that risk management is an inherently Bayesian endeavor, favoring notions of probability as a degree of belief given the knowledge at hand. The author argues that a “probability concept” is essential to describing and quantifying risks of any kind and gives an overview of how probability distributions can be used to describe uncertainty about a variable.
In addition to sources in the research literature, introductory texts such as Gelman et al. (10) and Christensen et al. (14) provide both practical and theoretical foundations for Bayesian statistics. In the following section, we offer Bayesian methods for estimation of resistance parameters DT, z, and Fo.
4. Estimation of Resistance Parameters
For estimation of a BI's DT value, the most common methods are the survival curve method and the fraction negative method. ISO-11138-1 recommends that both methods be used and that two DT values be reported. The Bayesian approach offers a framework by which the results from one method can be updated by the results of the second, such that one final estimate is provided for DT that incorporates the knowledge gained from both experiments. In this section, we consider estimation of DT using the survivor curve method, then we update this estimate using the Holcomb-Spearman-Karber (HSK) fraction negative procedure.
For the survival curve method, a population of microorganisms is exposed to successive short time periods of the sterilization conditions. ISO 11138-1 requires a minimum of five exposures. The surviving population at each exposure point is plotted in the log scale on the vertical axis versus time (in minutes) on the horizontal axis. From eq 3, it is clear that DT is simply the reciprocal of the slope of the survival curve, traditionally estimated using ordinary least squares (OLS) regression. If we let yi = logN(ti), β0 = logN(0), and 
, then motivated by eq 3, we write
where the errors, εi, i = 1, …, n, are taken to be independent and normally distributed with zero mean and variance σD2. From a Bayesian perspective, we require prior distributions for β0, β1 (or DT), and σD2. We treat these as independent and denote them by π(β0), π(β1), and π(σD2), respectively. We consider the construction of these priors in detail later. The joint posterior for the parameters of interest is obtained using Bayes' theorem:
where y and t denote the data vectors, l(·) is the likelihood function, and N denotes the normal distribution.
The HSK procedure for estimating DT was developed by Holcomb (15) as an extension of the Spearman-Karber methods used in analyzing serial dilution assays (16). The procedure is used to estimate UHSK, the mean expected time until a sample containing microorganisms becomes sterile from a set of experimental data in the quantal region (17). The quantal region is the exposure period over which a set of replicate test BIs exhibit a dichotomous response—some are positive for growth and the rest are negative (1). Assuming that the number of surviving organisms follows a Poisson distribution, the mean DT value, D̄T, is derived by Holcomb (15) as
where
Uk is the first exposure to show no growth of the replicates, d is the time interval between exposures, m is the number of replicates at each exposure, ri is the number of test samples showing no growth at exposure i, and k is the number of exposures. Note that eqs 9 and 10 are based on the limited version of the HSK procedure (6). If we let Xt denote the number of surviving microorganisms at time t, then Xt ∼ Poisson(N(t)) where N(t) is assumed to be the mean or expected number of microorganisms recovered after t min of exposure (see 15 and 18). It follows from eq 3 that N(t) = 10logN(0)−t/DT. Therefore, the probability of a sterile sample is
When m samples are tested for growth, the probability of obtaining w samples that are negative and (m − w) samples that are positive can be expressed using the binomial distribution as
From a Bayesian perspective, since we have posterior distributions for DT and N(0) from the survival curve method, we can simply use these distributions as the priors in the HSK procedure. That is,
where w denotes the data vector of negative responses, π(DT) and π(N(0)) are the lognormal approximations to the posterior distributions from the survival curve method, and Bin denotes the binomial distribution. Thus, we can obtain an up-to-date posterior distribution for the decimal reduction factor using the total data from two independent experiments.
From eq 3, it is evident that the value for the initial viable count, N(0), is a key factor in establishing the sterilization parameters of time and temperature that will result in an optimum process. However, enumeration techniques for determining the number of viable microorganisms are highly sophisticated and prone to measurement uncertainty. As such, it is appropriate to represent the uncertainty about N(0) with a probability distribution both before and after enumeration experiments. In order to ensure non-negativity, we can construct a mildly informative lognormal (LN) prior for N(0) where the median value is as close as possible to the reported viable count of a BI supplied by the manufacturer labeling, NM(0). The variance of the prior distribution can be chosen such that approximately 95% of the prior density is contained in the interval (0.5 × NM(0), 3 × NM(0)). This interval corresponds to the ISO-11138-1 allowance for the tested viable count of a BI to be between 50% and 300% of NM(0).
In terms of the decimal reduction factor, at a temperature of 121.1 °C, the typical D121.1 value for spore forming bacteria is in the range of 1 to 5 min (see 19–21). According to USP<1035> (22), Bacillus stearothermophilus, which is a commonly used BI for steam sterilization, has a typical D121.1 value of 1.9 min with a maximum of 3.0 min. In order to ensure non-negativity, we can construct a mildly informative lognormal prior for D121.1 where the median value is approximately 2.0 and the variance is chosen such that roughly 95% of the prior density is contained in the interval (0.5, 7.0). Finally, for σD, we assign a relatively non-informative U(0,A) prior where A is chosen sufficiently large as to prevent posterior influence.
The resistance coefficient, z, is traditionally estimated from multiple DT estimates obtained at different exposure temperatures. Given that a linear relationship between temperature and logDT is assumed (see Section 2), the estimator for z is the negative reciprocal of the OLS estimator for the slope coefficient of the regression equation relating T and logDT (23). From a Bayesian perspective, a posterior distribution for z can be directly induced using eq 5 once posterior distributions for D121.1 and DT are obtained.
Finally, the posterior distribution for Fo can be induced by summing over the n posterior distributions for the lethality rate in eq 7. Thus, a probability distribution representing everything we know about Fo is available for inference and decision making.
5. The Risk of Non-sterility
As described in Section 1, sterility assurance is defined as the probability of a single viable microorganism occurring on an item after sterilization. The PNSU measure is commonly used in the pharmaceutical industry to measure sterility assurance, and PNSU ≤10−6—or, equivalently, a 6 log reduction of the initial microbial population—is required to consider a system to be “sterile”. Using the Bayesian methods described in Section 3, the posterior risk of non-sterility can be directly quantified. That is, eq 4 can be used to induce posterior distributions not only for PNSU but also for the achieved log reduction. The resulting posterior distribution for PNSU can be used to quantify the probability or “risk” of non-sterility, denoted by α. If we let θ represent PNSU and λ represent the log reduction, then α can be calculated as follows:
or, equivalently,
where π(θ|data) and π(λ|data) denote the induced posterior densities for PNSU and the log reduction, respectively. Eqations 11 and 12 provide (1 − α) × 100% assurance that PNSU ≤ 10−6. In the following section, we present an example of empirical estimation of N(0), DT, z, and Fo values and the probabilistic conclusions that can be drawn regarding these parameters using the Bayesian approach. We also evaluate the posterior risk of non-sterility achieved by the sterilization process as defined by α.
6. Example
Suppose that DT values are being estimated using both the survivor curve and HSK methods. Table I provides the exposure times and population recovered, N(t), for two temperature points: 121.1 °C and 115 °C for the survivor curve method. The manufacturer's labeled value for the viable count of the BI is NM(0) = 2 × 106. Per Section 4, we construct the prior N(0) ∼ LN(μ = ln(2 × 106), ϕ = 0.5), where
Raw Data for DT Value Determination Using the Survival Curve Method
The resulting prior has a median value of 2,001,593 with 95% of its density contained in the interval (750,080, 5,342,276). The prior for D121.1 is LN(ln(2), 0.65) with a median value of 2.0 and 95% of its density contained in the interval (0.56, 7.14). For D115, we elicit expert opinion and consult the literature to construct a mildly informative LN(ln(6), 0.65) prior with a median value of 6.0 and 95% of its density contained in the interval (1.85, 19.38). Note that the constructed priors for D121.1 and D115 imply that an induced prior for z based on eq 5 can have intractable values where the denominator is zero. Although the induced posterior for z is also susceptible to this problem, we find that the experimental data update the posterior distributions for D121.1 and D115 sufficiently to prevent it (see Figure 2). Finally, for σD, we assign a relatively non-informative U(0, A = 2) prior. The chosen value for the constant A is made such that for A ≥ 2, the medians and 95% credible intervals for σD remain stable (refer to Figure 1).
Credible interval and posterior medians for σD at increasing values of A.
(A) Prior and posterior distributions for D121.1. (B) Prior and posterior distributions for D115.
Using OLS, the estimators for the decimal reduction times are D̂121.1 = − 1/β̂1 = 1.9 and D̂115 = 6.1. For the Bayesian approach, we obtain a posterior distribution for D121.1 and D115 based on the model in eq 7 using Markov-chain Monte Carlo (MCMC) simulation via the Stan package for the R software. A total of 50,000 iterations were used for inference after a burn-in of 10,000 iterations and a thinning factor of 2. Convergence was verified using traditional diagnostic methods as described in Gelman et al. (10). The posterior distributions for D121.1 and D115 are presented in Figure 2. The posterior means for D121.1 and D115 are 1.9 min and 6.1 min, respectively, with 95% credible intervals of (1.8, 2.2) and (4.7, 8.1).
The data for the HSK procedure is presented in Table II. A total of k = 7 exposures are performed with m = 20 samples per exposure. The time interval between exposures is fixed at d = 1 min for T = 121.1 °C and d = 3 min for T = 115 °C. Taking the posterior distributions for D121.1, D115, and N(0) from the survival curve method as priors, we execute MCMC simulation using 50,000 iterations for inference after a burn-in of 10,000 iterations and a thinning factor of 2. Posterior distributions for the parameters of interest are depicted in Figure 2, where significant updating of D121.1 and D115 is noted after observing the data from the HSK procedure. The posterior means for D121.1 and D115 are 2.0 min and 6.1 min, respectively, with 95% credible intervals of (1.9, 2.1) and (5.9, 6.3).
Raw Data for DT Determination Using the HSK Procedure
In terms of the resistance coefficient, z, the OLS estimator is ẑ = 12.1 °C. For the Bayesian approach, we use the posterior samples from D121.1 and D115 to induce a posterior distribution via eq 5. The resulting probability density is provided in Figure 3(a). The posterior mean for z is 12.7 °C with a 95% credible interval of (12.1, 13.4).
(A) Posterior distribution for z. (B) Validation cycle times and temperature.
Following the determination of a posterior distribution for z, suppose that the target Fo for the sterilization, process, F0(tgt), is 12 equivalent minutes. For simplicity, only exposure samples at temperature values greater than 100 °C are considered. The validation cycle is depicted in Figure 3(b), where the average exposure temperature is recorded after each minute for a total of 21 min (n = 21 and Δt = 1). Using the Bayesian method, we can induce posterior distributions for Fo real-time, as well as at the completion of the sterilization process. In Figure 4(a), we see that, over time, the probability distribution for Fo becomes more diffuse or less concentrated. This phenomenon makes sense, as we are accumulating “uncertainty” about Fo after each additional exposure sample. The posterior distribution for Fo upon completion of the cycle is depicted in Figure 4(b). The posterior mean is 13.3 equivalent minutes with a credible interval of (13.1, 13.4). Note that using the traditional method for calculating Fo results in a value of 13.1 equivalent minutes, which is close to the posterior mean. However, such an approach provides no assessment of uncertainty about the estimate.
(A) Real-time posterior distributions for Fo after each minute, from left to right. (B) Posterior distribution for Fo after completion of the sterilization process.
To calculate the risk of non-sterility, α, we first induce a posterior distribution for the log-reduction value via eq 4. The resulting density is depicted in Figure 5. Using the expression in eq 12, the risk of non-sterility is or approximately 0.5%. Thus, a posteriori, we can claim with 99.5% probability that PNSU ≤ 10−6. Such a probabilistic risk statement cannot be made using the traditional frequentist methods used to estimate PNSU and demonstrates the advantage of a fully Bayesian approach to the estimation of resistance parameters. That is, the sterilization process can be evaluated and validated based on an acceptance criterion requiring an upper limit on the risk of non-sterility. Such a limit for α can be chosen by the user based on the specific product or process under consideration as well as the clinical risks associated with non-sterility.
Posterior distribution for the log-reduction value achieved by the sterilization process.
7. Discussion
Development of sterilization processes is dependent on estimation methods that are sophisticated and prone to measurement uncertainty, particularly as they relate to enumeration of microorganisms. In recognition of the uncertainty associated with such measures, standards such as ISO-11138-1 allow for very large margins of error when testing is performed to verify their values. Rather than depending on such margins, it is possible to apply a fully probabilistic approach to estimation of sterilization process parameters using Bayesian methods. We have shown in this paper that by employing mildly informative priors, credible intervals for the parameters DT, z, and Fo can be obtained and used for inference as opposed to traditional point estimates. Rather than reporting multiple values for DT, the Bayesian approach permits straightforward updating of the posterior distribution based on experiments using the survivor curve and fraction negative methods. Note that combining the information gained from both experiments is not required for the Bayesian approach to work, and reporting two distinct estimates for DT using the two methods is still possible. In terms of the initial viable count of a BI, a probability distribution for N(0) seems an appropriate approach to incorporate the uncertainty about this quantity into sterility calculations. Verification testing performed by users of BIs can rely upon credible intervals to provide assurance of resistance characteristics and population levels. Through a simulated example, we have shown that a risk-based approach can be used to validate a sterilization process. By quantifying risk in terms of probability statements, a Bayesian approach allows manufacturers of sterile drug products to control and adjust the risk of non-sterility via posterior distributions. Of course, the BI approach is still necessary to validate the mathematical approach presented in this paper, and the latter should not be depended upon solely for verification of sterilization efficacy. Finally, it is important to note that, although this paper has focused a statistical approach to quantifying the risk of non-sterility, the ultimate aim of any risk-based method should be to evaluate the risk to patients, or the clinical safety. Indeed, an anonymous referee for this article noted that BIs are used because they are spore-forming bacteria with comparatively extreme resistance to moist heat. The mean resistance of most pathogenic bacteria is far less than that of BI species. Given this microbiological fact, real clinical risk is most accurately determined by considering the biologically derived lethality against the mean resistance of a panel of disease-producing bacteria. It is also important to consider the fact that sterilization cycles must preserve the usefulness of materials being sterilized, including those that are sterilized in their final container, which requires that the efficacy of the drug product be preserved.
Conflict of Interest Declaration
The author(s) declare that they have no competing interests.
Acknowledgements
The authors would like to thank the three anonymous referees for their valuable comments and instructive suggestions, which helped to greatly improve the article. At the time that this article was written, the corresponding author was employed by Allergan, Inc.
- © PDA, Inc. 2017
