INTRODUCTION
We aim to implement our NP trapping systems in different steps of the wastewater treatment process. There are several factors that will affect the NP trapping efficiency for the proteorhodopsin (PR) bacteria and biofilm models. Our PR bacteria would be added to aeration tanks, where water movement is fast and turbulent, while our biofilm (attached to biocarriers), would be placed in the clarifier or sedimentation tanks, where water movement is calmer to keep the biofilm intact.
Due to the lack of literature on our proposed NP-trapping techniques using PR and biofilm, experimental trials and our prototype design were integral to the modeling process. Experimental trapping rates from our prototype were used to fit our model to the current trapping abilities of our PR construct (BBa_K2229400) and our biofilm construct (BBa_K2229300). After experimentally determining the rate constants for our PR bacteria and biofilm constructs, the mathematical models can be used to determine two objectives, given an initial NP concentration and a final target NP concentration.
Objective 1: What PR bacteria concentration is needed in the aeration tanks?
Objective 2: How many biofilm-coated biocarriers are needed in the sedimentation tank?
PROTEORHODOPSIN TRAPPING MODEL
Proteorhodopsin and citrate binding modeled as a ligand-receptor interaction
To model the binding of PR bacteria to CC-NPs, we used a coarse-grained model for ligand-receptor interaction (Ruiz-Herrero et. al 2013). The model is based on the chemical interaction between a freely diffusing ligand L (nanoparticle), and a cell membrane receptor R (PR bacteria), which combine to form a complex C (nanoparticle-loaded PR bacteria) in the following reaction scheme(Figure 4-1):
Figure 4-1 Reaction scheme for PR bacteria and nanoparticles. Our bacteria (R) and nanoparticles (L) bind with the affinity rate k(on) to form the complex C. Conversely, starting with the complex C, nanoparticles fall off bacteria with the dissociation rate k(off). Figure: Justin Y.
where kon is the binding rate constant of our PR bacteria to CC-NPs and koff is the rate constant of NPs dissociating from PR bacteria. L, R, and C are all functions of time because our PR bacteria bind to NPs over time, which decreases the concentration of free NPs (L) and available PR bacteria (R) while increasing the concentration of NP-loaded PR bacteria (C). Therefore, we can use the following differential equation to model the progression of NP trapping over time by our PR bacteria:
Equation 1: Binding and Dissociation Model
To determine the kon and koff rate constants, we fitted our model to match experiments (figures 4-3 and 4-4) that show a decrease in CC-AgNPs when mixed with a known concentration of our PR bacteria (BBa_K2229400). Determining kon and koff rate constants will enable us to use equation 1 and inform WWTPs what concentration of PR bacteria is needed to treat their tanks given a starting NP concentration and desired final NP concentration. Below, we will explain how we obtained both the kon and koff rate constants.
Determining NP binding rate (kon) and dissociation rate (koff) constants using experimental data
Binding-Only Model (kon only)
We first determined kon while assuming a best-case scenario where NPs do not fall off of PR once they bind, which means that koff is zero. Thus, our initial model can be described with the following equation:
Equation 2: Binding-Only Model
Since the time-dependent functions [L] and [R] for binding between CC-NPs and PR bacteria are unknown in literature, we used finite-difference methods (FDM) to model this equation. (Click here to learn more about how we used FDM!)
We obtained kon, the binding rate constant, from experimental data where CC-AgNPs were mixed with a known concentration of our PR bacteria (BBa_K2229400) (figures 4-2 and 4-3 below).
Figure 4-2 The binding rate kon for our PR bacteria was determined experimentally. A known concentration of our PR bacteria (BBa_K2229400) was mixed with a known concentration of citrate-capped AgNPs, and the decrease in nanoparticle concentration over 5 hours was measured. By inputting the known bacteria and nanoparticle concentrations, we fitted our model to our experimental data to find kon. Figure: Justin Y.
Figure 4-3 Determining kon. We obtained kon, the binding rate constant, from experimental data where CC-AgNPs were mixed with a known concentration of our PR bacteria (BBa_K2229400). Experiment: Justin Y
To find the kon value for our PR bacteria, we fitted our model to our experimental data (Click here to learn how we fit our model using RMSE!). 1.9×10-7 µL cells-1 hr-1 was the kon value for which the Root Mean Squared Error (RMSE), a measure of the error between the model and the experiment, was the lowest. However, as Figure 4-3 shows, the Binding-Only Model does not reflect how the experimental data reaches an asymptote towards the end of 5 hours. One of the things that may account for this discrepancy between our Binding-Only Model and our experimental data is NPs falling off of the PR bacteria. As we show below, our revised model does a better job of describing the experimental data after accounting for this effect. Since our Binding-Only Model does not take the NP dissociation rate constant, koff, into account, we now needed to determine this value.
Binding and Dissociation Model (kon and koff)
For the first hour, few NPs are bound by the bacteria and there are many PR receptor sites open, so we assumed that the number of NPs falling off is negligible in the first hour. Thus, we assume that the dissociation rate, koff, is 0 from t = 0 to t = 1. Based on this assumption, we found a kon value that fitted our model to the first two points (the first hour) of our experimental data. As the yellow curve (step 2, figure 4-4 below) shows, the model fits the first two data points but falls below the experimental data as time goes on. We assumed that this discrepancy is due to NPs falling off the bacteria. Thus, our next step was finding a koff value that would fit our model to the asymptotic behavior of our experimental data. 0.32 hr-1 was the koff value that yielded the lowest RMSE. (Click here to learn how we found koff using RMSE!). The green curve (step 3 figure 4-4 below) includes this koff value, and fits our experimental data much better with the lowest RMSE. Thus, our final Binding and Dissociation model includes both kon and koff.
Figure 4-4 Determining Koff to better fit our model to the experimental data. Binding-only model (blue) does not reflect how the experimental data reaches an asymptote at the end of 5 hours. We assumed that this discrepancy is due to NPs falling off of bacteria. To account for this, we determined a new Kon that fits the first two points of experimental data (yellow), and then we fitted the rest of the curve using a dissociation constant (Koff). The resulting curve (green) including Kon and Koff fits our experimental data much better.
CC-NP Trapping by PR Bacteria Calculator
We developed two calculators to help WWTPs use our PR bacteria to clean up CC-NPs.
Calculator 1 allows WWTPs to input their NP concentration, their target NP concentration, and the time water spends in the tank to determine the initial PR bacteria concentration they need to add.
Calculator 2 allows WWTPs to input their NP concentration, PR bacteria concentration they plan to add, and the time water spends in the tank to determine the final NP concentration of the water leaving the tank.
Example Application of Completed Model
We used our Binding and Dissociation model (figure 4-5) to determine the trapping of CC-NP concentration over time by our PR bacteria using the kon and koff values above (3.5×10-7 µL cells-1 hr-1 and 0.32 hr-1 respectively) and time intervals of 0.1 hours (click here to learn why).
Figure 4-5 Applying our model to our experimental conditions.
We used our Binding and Dissociation model to predict NP trapping for initial conditions that matched those of our experimental trial. The final NP concentration that our model predicts was 0.693 µM, which has a 2.19% error compared to our experimental value, 0.708µM.
In this example, the initial conditions of L and R were set to the same values as our experimental trial, which means that [NP] = 1.078 µM and [PR bacteria] (at OD600=0.712) = 569600 cells/µL. Under these conditions, our model predicts that NP concentration after 5 hours is 0.693 µM (the percent difference of our modeled value from our experimental value (0.708 µM) is 2.19%). WWTPs can obtain a graph like this one by inputting the variables specific to their treatment plant, such as initial NP concentration and how much time water spends in the tank.
Calculation Explanations for FDM and RMSE
Finite-Difference Method (FDM) Explanation
We used Euler’s Method to create our FDM model, and our formulas for the Binding-Only Model are shown below (table 4-1). As the table shows, by plugging in values for L and R at t = 0, a fixed value for kon, and the time interval Δt, dC/dt = Y1 at t = 0 can be calculated. In other words, we can find the trapping rate at t = 0. As per FDM methods, we assume that the trapping rate is Y1 for t = [0, t0 + Δt] = [0, t1], so 
= [NP] trapped during the time interval t = [0, t1]. Therefore, to find the new L at t = t1, we subtract (Δt)Y1 from L0 to get L1. The same process is done for R to find R1 because the decrease in [NP] is directly proportional to the decrease in unsaturated bacteria. Then, using the new values L1 and R1 and the same values of kon and Δt, a new trapping rate Y2 for t = [0, t1 + Δt] = [0, t2] is calculated. This process is repeated for each subsequent value of t.
Table 4-1 Binding-Only Model Formulas These are the formulas we entered into Excel to run Equation 2 of our model (where koff = 0) using FDM.
Table 4-2 shows a sample of our calculations for the Binding-Only Model using the formulas above, initial L and R conditions that match our experimental trial, the value of kon determined from experimental data, and Δt = 0.1 up to 1 hour.
Table 4-2 Binding-Only Model Calculation Sample
This table shows the values we got for initial conditions matching our experimental trial, our value of kon, and time intervals of 0.1 hours going up to t = 1 hour.
At t = 0, no NPs have been trapped and all PR bacteria are available, so we can set L = initial [NP] and R = initial [bacteria], and use these initial conditions to calculate subsequent values.
For our Binding and Dissociation Model, the formulas were modified to include the effect of koff, which is shown below (table 4-3). This shows one of the advantages of Euler’s Method and other FDM models ‒ the form of dC/dt is adaptable to further variables or complications, and the model can still provide results for complex equations.
Table 4-3 Binding and Dissociation Model Formulas
These are the formulas we entered into Excel to run Equation 1 of our model (which incorporates both kon and koff) using FDM.
Table 4-4 shows a sample of our calculations for the binding and dissociation model using the formulas above, initial L and R conditions that match our experimental trial, the values of kon and koff determined from fitting to experimental data, and Δt = 0.1 up to 1 hour.
Table 4-4 Binding and Dissociation Model Calculation Sample
This table shows the values we got for initial conditions matching our experimental trial, our value of kon and koff, and time intervals of 0.1 hours going up to t = 1 hour.
When using FDM, it is important to keep in mind that values after the set initial conditions are numerical approximations. This means that the larger the value of Δt, the more inaccurate the final result, as shown in figure 4-6.
Figure 4-6 Effect of Changing Δt on FDM Calculations
Decreasing Δt means more time interval divisions for a given time frame. This results in a smoother and more accurate curve. Thus, the smaller Δt is, the more accurate the model.Figure: Justin Y.
For example, Figure 4-7 below shows how the trapping curve changes with different values of Δt. Larger values of Δt result in a final [NP] that is too low. This is because trapping rate is fastest at the beginning (when there are more unbound NPs and unsaturated bacteria), and large values of Δt assume this faster trapping rate for a larger interval of time, thus resulting in a faster decrease in [NP] than if smaller values of Δt were used. To account for error caused by FDM, we ran our model with decreasing values of Δt until there was a negligible difference in final [NP], which was at Δt = 0.1. Since the error due to FDM is negligible for Δt = 0.1, our model used time intervals of Δt = 0.1 hours.
Figure 4-7 Determining an appropriate Δt for our model
We tried four values of Δt (1, 0.5, 0.25, and 0.1), and each resulting curve is shown in the figure. The red curve, which has the largest time intervals (Δt = 1) is the least smooth while the blue curve, which has the smallest time interval (Δt = 0.1) is the smoothest. There is a negligible difference between the green and blue curve on the graph, which means that we can use Δt = 0.1 for our model and assume that error caused by FDM is negligible.
Minimizing Root Mean Square Error (RMSE) to find kon and koff
To determine the best fit of our model to our experimental data, we calculated the Root Mean Square Error (RMSE). We did this by squaring the difference between each experimental value and its respective modeling value, averaging these values, and finding the square root. The lower this RMSE value, the better the fit of our model to our experimental data.
For our binding-only model (the blue curve in Figure 4-4), the RMSE was calculated for all of the experimental data points, and we found an initial kon value that gave the lowest RMSE. For our binding and dissociation model, we first found a new kon value that fit the first two experimental data points, then found a koff value that fit the rest of the data points. To find the kon value, the RMSE was calculated for the first two points. To find the koff value, the RMSE was calculated for all the points.
BIOFILM TRAPPING MODEL
Free-floating PR bacteria and biofilm have different characteristics, so we have to consider different variables in our biofilm trapping model. Our goal here is to use mathematical modeling to 1) understand how biofilms interact with NPs, and 2) use this information to predict how many biofilm-coated biocarriers should be used in WWTP sedimentation tanks.
Based on experiments described below, we first determined that biofilm surface area (SA) is the main factor affecting NP trapping rate. This led us to develop a biofilm trapping model which treats NPs as solutes in wastewater that get trapped as they are carried into the biofilm surface by the movement of water (as shown in Figure 4-8).
Figure 4-8 Variables affecting how biofilms trap NPs NPs in wastewater get trapped as they contact a biofilm surface. Therefore, NP trapping depends on the speed of water carrying NPs into a biofilm's surface and the surface area of biofilm. Figure: Justin Y.
We modeled this interaction using a differential equation (click here to see how we derived this equation):
which relates the change in NP concentration to an effective trapping rate constant (ktrap), the speed of the water carrying the NP into the biofilm (Swater), the total surface area of biofilm (SAbiofilm), and the total volume of the treatment container (Vwater total).
As with the PR bacteria model, we used FDM to solve this equation and ran experiments to find ktrap based on experimentally controlled values of the other variables.
Finally, we developed two calculators to help WWTPs use our biofilm to clean up NPs.
Determining Key Parameters of Biofilm-NP Interaction
To model our biofilm’s NP trapping rate, we first needed to determine the significance of factors that would affect the biofilm-NP interaction. Initially, we assumed that the surface area (SA) and volume of the biofilm would be the two major factors that affect trapping rate since nanoparticle trapping depends on physical interaction between biofilm and nanoparticles.
Biofilm Volume vs NP Trapping
In our experiments, we observed that NPs were trapped in a matrix of biofilm (see our experiment), but we did not know how deep the NPs would travel into the biofilm layer. When we tested 4 different volumes of biofilm for their NP trapping abilities over a period of 5 hours (to mimic WWTP sedimentation tank conditions), we found that volumes greater than 1-1.5 mL—which covered the entire bottom of the cylindrical container—did not significantly improve NP trapping (see our experiment). This led us to conclude that NPs are not likely to travel through the entire biofilm layer.
Biofilm Depth vs NP Trapping
Wanner et al. (2006) previously reported that biofilms could only trap particles up to a depth of 200 μm. To see if our own biofilm constructs have a similar trapping depth, we analyzed our volume trial results by converting the biofilm volumes to depths. Since our biofilm is different from the ones used by Wanner et al., it might have a different maximum depth for trapping NPs.
In our experiments, biofilm was added directly to the cylindrical containers, where volume can be found by the equation v=πr2h. We tested 4 different biofilm volumes, so the height (h) or depth of the biofilm was found by plugging in known values for volume (v) and radius (r).
Figure 4-9 Effects of biofilm depth on NP trapping When biofilm depth exceeds a certain limit, depth no longer impacts NP trapping. We experimentally determined that a depth exceeding 1000 μm does not increase our biofilm's NP trapping ability. Thus, we made sure to use a biofilm depth over 1000 μm in our future experiments.
We found that more AuNPs were trapped when the biofilm increased from 500 to 1000 μm, but when the depth exceeded 1000 μm, the amount of AuNPs trapped by our biofilm remained similar (figure 4-9). Therefore, our experimental results agree with the general conclusion of Wanner et al. (2006) that when biofilm depth exceeds a certain limit, depth no longer impacts NP trapping. In future experiments, we made sure to use biofilm with a thickness over 1000 μm (1 mm) to ensure that the depth is not a factor affecting NP trapping.
Biofilm Surface Area vs NP Trapping
Our biofilm volume and depth experiments allowed us to conclude that, once biofilm depth exceeds 1000 μm, any additional depth does not affect NP trapping, suggesting that surface area might have a larger effect on trapping rate. To test this, we ran an experiment using two different containers. We assumed that the biofilm surface area exposed to the AuNP solution is equal to the surface area of the containers’ bottoms (10.18 cm2 and 2.01 cm2). In both containers, biofilm depth exceeded 1000 μm. We found that the higher biofilm surface area resulted in a higher rate of AuNP trapping (figure 4-10).
Figure 4-10 Increased biofilm surface area leads to a higher rate of AuNP trapping We ran two identical experiments using different biofilm surface areas. Our results show that larger surface areas increase AuNP trapping.
Since a larger surface area will trap NPs at a faster rate, we designed biocarriers with a large surface area. Knowing the biofilm surface area needed and the surface area of our biocarrier, we can calculate the number of biocarriers needed to trap a known amount of NPs.
Modeling Biofilm’s Ability to Trap NPs
As long as biofilm depth exceeds 1000 μm, surface area is the main factor that affects NP trapping. We model how surface area affects the NP trapping rate when NPs come into contact with our biofilm. Our model treats NPs as well-mixed solutes that are carried into the biofilm by water. Due to the small mass of NPs, we assume that the speed of NP movement is the same speed as the turbulent water. Hence, knowing the volume of water in contact with our biofilm (Vwater) and NP concentration in the water ([NP]), we can determine how many NPs are in contact with our biofilm.
Since Vwater is the product of biofilm surface area (SAbiofilm), time elapsed (T), and the speed of water in contact with biofilm(Swater), our equation becomes:
To model how many NPs in contact with biofilm are actually trapped, we added a trapping rate constant ktrap, which is a value between 0 (none trapped) and 1 (all NPs in contact with biofilm are trapped):
In our model, we also assume that NPs trapped by biofilm would sink into the biofilm matrix, and remain embedded. Thus, we do not have an expression for dissociation of trapped NPs from the biofilm.
Taking the volume of a sedimentation tank (Vtank) into account, we can convert NPtrapped in biofilm into NP concentration of water in the sedimentation tank([NP]). The differential equation below models the change in NP concentration in the tank over time:
We determined ktrap based on experiments with known values for Swater, SAbiofilm, [NP]initial, and Vtank. In these experiments, AuNP solution and a layer of biofilm was added in plastic cylinders to simulate a WWTP sedimentation tank. With a central rotor slowly rotating, we measured the absorbance of the AuNP solution every 20 minutes for a total of 4 hours (see Modeling in our lab notebook). We then fit our model to our experimental data to find ktrap, with the same RMSE method used in our PR model (Figure 4-11; link to biofilm’s RMSE method).
Figure 4-11 Determining ktrap, the trapping efficiency rate constant. To fit experimental data where AuNP solution interacted with our biofilm, we used RMSE method to minimize error.
Biofilm-NP Trapping Calculator
Using our model, we developed two calculators to help WWTPs use our biofilm to clean up NPs.
Calculator 1 allows WWTPs to input their initial NP concentration ([NP]i), targeted final NP concentration ([NP]f), speed of water in contact with biofilm (Swater), total volume of the tank(Vtank), and the total time water spends in the tank (T) to determine the surface area of biofilm (SAbiofilm), which we convert to number of biocarriers needed to reach the targeted final NP concentration.
Calculator 2 allows WWTPs to input their initial NP concentration ([NP]i), speed of water in contact with biofilm (Swater), total volume of water (Vtank), time water has to spend in the tank (T), and # of biocarriers to determine the final NP concentration.
Example Application of Completed Model
In figure 4-12, we used our model (with a time interval of 5 min-- click here to learn why) to determine the decrease in [NP] in the water over time using the ktrap rate constant of our biofilm. In this example, the initial conditions were set to the same value as our experimental trial:
- [NP]i of 24.628 µM
- Swater of 66.72 cm/min
- Vtank of 0.310 liters
- T of 24 hours
- SAbiofilm of 29.22 cm2
Figure 4-12 Change of NP concentration over time Using our biofilm trapping differential equation, we modeled the change of NP concentration over 24 hours using our biofilm.
Under these conditions, our model predicts that NP concentration after 24 hours is 5.104 µM (the percent decrease in NP concentration is 79.3%).
Calculation Explanation for FDM and RMSE
Finite-Difference Method (FDM) Explanation
As with the PR model, we used Euler’s method to approximate our differential equation. (link to Euler’s method) The main difference between our biofilm FDM model and the PR FDM model are the variables used. For our biofilm FDM model, we plugged in the variables: the speed of water Swater, total volume of water in the tank, the biofilm surface area, a fixed ktrap and the initial NP concentration at time t=0. We also used a shorter time interval of Δt = 5 minutes to get a higher accuracy (Figure 4-6). As per FDM methods, we assume that the rate of decrease for NP concentration is approximately constant in each 5 minute interval. The new rate of NP concentration change is determined after each interval according to our biofilm-NP trapping differential equation using the most recent NP concentration.
Table 3-5 shows a sample of our calculations for the biofilm-NP trapping model in the first hour using initial conditions that match our experimental test. The value of ktrapused here is the most accurate trapping efficiency rate constant as found by our RMSE method (link to RMSE).
Minimizing Root Mean Square Error (RMSE) to find ktrap
To determine the best fit of our model to our experimental data, we calculated the Root Mean Square Error (RMSE). We did this by squaring the difference between each experimental value and its respective modeling value, averaging these values, and finding the square root. The lower this RMSE value, the better the fit of our model to our experimental data. For our biofilm-NP trapping model, the RMSE was calculated for all of the experimental data points, and we found a ktrap value that gave the lowest RMSE.
Experimentally Measuring Swater
Dye-Drop Experiment Explanation
Our biofilm-NP trapping model depends on the speed of water in contact with our biofilm Swater. In our experiment, the speed of the water is due to mixing by the central rotor. Since we do not have a precise sensor that can measure this, we designed an experiment using food coloring to determine the speed of water in this setup. We placed biofilm at the bottom of the container, so Swater should be equal to the vertical speed of water in the tank.
We carefully dropped food coloring from the surface of the water, and measured the time it took for any part of the dye to reach the bottom of the container. Dividing the vertical distance of the tank by this elapsed time gave us an approximate Swater.
REFERENCES
Ruiz-Herrero T, Estrada J, Guantes R & Miguez DG (2013). A Tunable Coarse-Grained Model for Ligand-Receptor Interaction. PLoS Comput Biol 9(11): e1003274. doi:10.1371/journal.pcbi.1003274
Wanner, O., Ebert, H.J., Morgenroth, E., Noguera, D., Picioreanu, C., Rittmann, B.E., Van Loosdrecht, M.C.M., 2006. Mathematical modeling of biofilms. IWA Scientific and Technical Report No. 18, IWA Task Group on Biofilm Modeling
