do not necessarily reflect the views of UKDiss.com.
Evaluation and Mitigation of the Temporal Evolution of Microbial Contamination Risk in Surface Water Systems
Abstract
This thesis presents a methodology for evaluating the temporal evolution of microbial contamination risk in surface water systems. The objective of this research is to advance the field of risk assessment across the realm of environmental engineering by creating a methodology that is applicable to numerous complex environmental phenomena, as well as developing the tools and theories necessary to make this methodology pragmatically viable. This is accomplished in part by developing a probabilistic mass conservation model to describe large scale dynamics of microbial contaminants from nonpoint sources throughout a watershed, in a manner that embraces the uncertainty inherent to these phenomena. With a developed microbial fate and transport model that attributes contamination to discrete spatial sources and pathways, the action that must be taken to minimize this contamination risk can be determined. Any method utilized to minimize nonpoint source microbial contamination risk will have an associated cost of implementation, so an optimization algorithm is applied to the results of the probabilistic mass conservation model to determine the optimal placement of control measures throughout the watershed necessary to achieve an objective reduction in contamination risk. To aid in the practicality of risk assessment within complex environmental systems, a rapidly converging algorithm is developed to efficiently propagate uncertainty through complex systems, such as that which is described by the probabilistic mass conservation model. This describes the time and space evolving probability density function of microbial contamination, which in turn is used to determine how risk changes throughout time. Finally, a standardized methodology for parameterizing the probabilistic mass conservation model is presented, which in turn allows for the derivation of certain transition rates using statistical and physical first principles. The ability of environmental engineers and watershed managers to perform a large-scale, systems-driven assessment of risk in environmental systems has not been previously viable, due to the irreducible complexity of environmental processes and the significant computational cost demanded by conventional means of uncertainty propagation for risk assessment. This project provides a versatile, scalable platform for the assessment of the temporal evolution of risk in environmental systems, applied specifically to microbial contamination in surface water systems.
Preface
In Chapter 1, a methodology is presented to determine the temporal evolution of microbial contamination risk in surface water systems. To describe the temporal evolution of risk of microbial contamination, system components must first be identified. System components are the potential sources and pathways of microbial contamination that contribute to water quality failure. These components are identified by first relating each spatial segment of the watershed to its contribution of microbial contaminants to the surface water system using sets of directed paths. This captures the influence that spatially dependent characteristics of the watershed have on nonpoint source microbial contamination. Next, the stochastic spatial-temporal behavior of microbial contamination is described within each component using a probabilistic mass conservation model. Finally, the system components must be recombined within the system reliability framework to create the metric of temporal evolution of risk. Using the temporal evolution of risk, sustainability that aid in the interpretation of the temporal evolution of risk metrics will be created.
In Chapter 2, a linear optimization algorithm is applied to results of probabilistic mass conservation model to determine the optimal implementation of control measures and the reduction of contamination risk. This can be summarized in three different schemas: single point control, multiple point control, and multiple point control through time.
In Chapter 3, an algorithm is developed to efficiently propagate uncertainty through complex systems, such as that which is described by the probabilistic mass conservation model. The algorithm allows the description of time and space evolving probability density function of microbial contamination, which in turn is used to determine how risk changes throughout time. The general function of this algorithm is to determine an output uncertainty by identifying the combination of input values that result in an output that is one standard deviation above the mean output, where the output standard deviation is unknown. This is accomplished by numerically approximating the influence that each parameter has on the output, then using this influence in tandem with a converging algorithm to determine the standard deviation of the output.
In Chapter 4, a methodology is presented for parameterizing the probabilistic mass conservation model, specifically through deriving transition rates using statistical and physical first principles. The importance of treating transition rates in this manner is demonstrated by showing that a number of transition rates can be derived using statistical first principles to determine the conditional transition probability. It is shown that while such derivations provide elegant formulas for certain transition rates, they can become extremely complex for other parameters. However, this increase in complexity does not preclude the importance of treating transition rates as conditional probabilities divided by a time increment. The Bayesian nature of these transition rates has a bearing on the empirical approximation of transition rates that have unwieldly purely theoretical derivations, and provides insight into the probabilistic distribution of these transition rates.
This page Left Intentionally Blank
Chapter 1: Temporal Evolution of Risk and Sustainability Metrics
1.2.1 Reachability and Directed Paths
1.2.2 Microbial Transport along Directed Paths
1.2.4 Temporal Evolution of Risk
1.2.5 Sustainability as a Function of Time
1.3.1 Temporal Evolution of Risk Results
1.3.4 Sustainability Index Results
1.4.1 Combination of Metrics for Engineering Judgement
1.4.2 Use of Metrics in Decision Making
1.4.3 Potential Applications in System Controllability
Chapter 2: Contributor Based Control of Microbial Contamination
2.2 Proportional Reduction in Contributor Based Control
2.2.1 Proportional Reduction in Single Point Control Schema
2.2.2 Single Point Control Application
2.2.3 Multiple Point Control Schema
2.2.4 Multiple Point Control Application
2.2.5 Multiple Point Control through Time Schema
2.3 Control Application Summary
2.3.2 Possible Means of Control through Application of Management Practices
2.3.3 Grass Buffer Strips Effectiveness
2.3.4 Grass Buffer Strips Feasibility
2.3.5 Source Control Effectiveness
2.4.1 Single Point Control Results
2.4.2 Single Point Control Discussion
2.4.3 Multiple Point Control Results
2.4.4 Multiple Point Control Discussion
2.4.5 Multiple Point Control through Time Results
2.4.6 Multiple Point Control through Time Discussion
Chapter 3: Rapidly Converging Uncertainty Propagation Algorithm
3.2 Uncertainty Propagation Theory
3.4 Uncertainty Propagation with Correlated Random Variable Inputs
Chapter 4: Parameterization of Probabilistic Mass Conservation Model
4.2 Probabilistic Nature of Transition Rates
4.3 Derivation of Transition Rates
4.3.1 Lateral Transport Transition Rate Substantiation
4.3.2 Settling Transition Rate Substantiation
4.3.3 Uptake Transition Rate Substantiation
4.3.4 Empirical Derivation of Parameters Methodology
4.4 Probabilistic Distribution of Transition Rates
Appendix A: Probabilistic mass conservation model MATLAB Code
Appendix B: Hill Slope and Stream Tables
Appendix C: Temporal Evolution of Risk MATLAB Code
Appendix D: Control Optimization MATLAB Code
Appendix E: Rapid Uncertainty Propagation Algorithm MATLAB Code
Appendix F: Temporal Sustainability Index MATLAB Code
Index of Figures
Figure 1.1.A: Example of Simple Watershed
Figure 1.1.B: Extracted Dendritic Network
Figure 1.2: Schematic of Resulting System Components
Figure 1.3: Microbial Transport Model
Figure 1.4: Microbial Transport Finite Volume Method
Figure 1.5: Fault Tree Diagram
Figure 1.6: Collapsed Fault Tree Diagram
Figure 1.7.A: Shepherd creek watershed
Figure 1.7.B: Sub-Catchment of Shepherd Creek Watershed
Figure 1.8.A: Concentration Mean Curve of System Components
Figure 1.8.B: Concentration Standard Deviation Curve of System Components
Figure 1.9: Concentration Mean Curve of System with Measured Data
Figure 1.11: Conventional Resilience
Figure 1.13: Vulnerability Results
Figure 1.14: Sustainability Index Results
Figure 2.1: Sample Watershed with Objective Points
Figure 2.2: Concentration Mean Curve at Point B
Figure 2.3: Concentration Mean Curve at Point C
Figure 3.1: Convergence of Guessed and Actual Standard Deviation
Figure 3.2: Convergence of Uncertainty Propagation Algorithm
Figure 3.3: Comparison of Uncertainty Propagation Algorithm and MCS for Time Dependent System
Figure 3.4: Mean Solution of Lorenz Equation
Figure 3.5: Comparison of Uncertainty Propagation Algorithm and MCS for Time Dependent System
Figure 4.1: Free Body Diagram of Soil Particle on Hill Slope
Figure 4.2.A: Sampling Results from Beta Distribution with
clear x
for it=1:length(std_array_new)
%Evaluate the Mean and standard deviation from each component at each point in time
mean_at_time=mean_array_new(it,:);
std_at_time=std_array_new(it,:);
%Evaluate the combined mean and standard deviation from all components at
%each point in time
m=sum(mean_at_time);
s=sqrt(sum(std_at_time.^2));
%This improves stability for very low probability of failures
if s<1e-11
Prob(it)=0;
else
%Determine parameters of lognormal distribution
mu = log((m^2)/sqrt(s^2+m^2));
sigma = sqrt(log(s^2/(m^2)+1));
%Determine probability of failure through cumulative lognormal
%distribution.
Prob(it)=1-logncdf(Limit,mu,sigma);
end
end
Appendix D: Control Optimization MATLAB Code
%Controlability
clear all
clc
N=11; %Number of Nodes
C=9; %Number of Control Points
%Import Contamination Data for Control Points based on Probabilistic Mass
%Conservation Model
L=transpose(xlsread(‘Shep Control.xlsx’,’Sheet3′,’B8:J18′));
%This is our goal load at the point of interest, so the accumulation of all
%the individual nodes should be below this.
Lg=5e4;
%Proportional reduction
%So this portion creates the matrix we are going to use. The decision
%variable X is the proportion that we reduce each node load by.
for ic=1:C
for i=1:N
A(ic,i)=-L(ic,i);
end
end
%Establish Initial Constraint Matricies
AB=zeros([2*N+C N]);
AB(C+1:N+C,:)=eye(N,N);
AB(N+C+1:2*N+C,:)=-eye(N,N);
%Establish constraint that sum of the loads minus their reduction to be less than the
%goal.
for ic=1:C
AB(ic,1:N)=A(ic,1:N);
B(ic)=Lg-sum(L(ic,:));
end
%Determine Maximum proportional reduction
B(C+1:N+C)=0.65;
B(N+C+1:2*N+C)=0;
B=transpose(B);
z(1:N)=1;
%Perfom Linear Optimization to determine optimal proportional reduction
[X Z]=linprog(z,AB,B)
Appendix E: Rapid Uncertainty Propagation Algorithm MATLAB Code
%Rapid Uncertainty Propagation Algorithm
clear all
%Establish Inputs
s=rng(3);
Vars=4;
InputMean=rand(Vars,1)*1;
Inputstd=rand(Vars,1)*0.1;
%Initialize count
Count=1;
%Specify Time of Evalutation (if relevant)
eval_time=500;
%Create Vector of Input Means
In=InputMean;
%Run function with output gC
GcFun
Count=Count+1;
%Because the mean inputs were entered into the model, the output is the
%mean output
Mean=gC;
%Specify Variable Number
var_num=length(InputMean);
%Create range of input weights from -3 to 3, create matrix were each weight
%is varied across this range while all others are 0.
gindex=1;
range=3;
clear w
for q=1:var_num
for iw=-range:range
w(gindex,1:2)=zeros(1,2);
w(gindex,q)=iw;
gindex=gindex+1;
end
end
%Sweep across the input weight matrix
for iW=1:length(w)
%Apply input weights to inputs
In=InputMean+Inputstd.*transpose(w(iW,:));
%Run Model
GcFun
Count=Count+1;
%Save Results
DiffOutputs(iW)=gC;
InSave1(:,iW)=In;
end
y=DiffOutputs;
%Find vector of differences between output and mean
D=(DiffOutputs-Mean);
%For stability
for itna=1:length(D)
if isnan(D(itna))==1
D(itna)=0;
end
end
%Determine Response Space Characteristics
syms x
iRc=1;
%Determine Dimension Size
Dimsize=length(w)-6;
for iw=1:7:Dimsize
%Perform linear regression on response space
[polyvec]=polyfit(InSave1(iRc,iw:iw+6)-InputMean(iRc),D(iw:iw+6),1);
%Save first parameter of linear regression
Mvec(iRc)=polyvec(1);
%Determine correlation coefficient between realizations across
%characteristic of response space
M=corrcoef(-3:1:3,D(iw:iw+6));
%Compare to correlation coefficient between differences across
%characteristic of response space
MP=corrcoef(InSave1(iRc,iw:iw+6)-InputMean(iRc),D(iw:iw+6));
%Extract Correclation Coefficient
Rc(iRc)=M(2,1);
RcP(iRc)=MP(2,1);
%Determine the influence vector
Infl(iRc)=mean(diff(D(iw:iw+6)));
%Indexing
iRc=iRc+1;
end
%Allow for iteration of guesses
for iCC=1:1
%Save initial Influence guess
Infl_f=Infl;
%Solve for adjusted influence
iW=1;
%Clear as necessary
clear w
clear A DiffOutputs StdGuess
DiffOutputs(1)=gC;
%Create initial Guess
StdGuess(1)=Mean;
%Save Guess
Std_Record(1)=StdGuess(1);
Esave=1;
%Determine convergence threshold
Thresh=Mean*1e-04;
%While system is above this threshold…
while abs(DiffOutputs(iW)-Mean-Std_Record(iW))>Thresh
%Based on guessed influence and guessed STD, generate vector of guessed
%input weights
w=Infl./StdGuess(iW);
%Apply guessed input weight vector to inputs
In=InputMean+Inputstd.*(w(:));
%Run Model
GcFun
%Save Output
Count=Count+1;
DiffOutputs(iW+1)=gC;
%Evaluate Local response space gradients about the guessed input weights
iRc=1;
for iw=1:7:Dimsize
%Spline interpolate locally
NewD(iw:iw+6)=interp1(-3:3,D(iw:iw+6),[w(iRc)+3*w(iRc)/-10]:w(iRc)/10:w(iRc)+3*w(iRc)/10,’spline’);
Infls(iRc)=mean(diff(NewD(iw:iw+6)))/w(iRc)*10; %Local RS gradient
iRc=iRc+1;
end
%Record Outputs
Std_Record(iW+1)=StdGuess(iW);
w_record(:,iW+1)=w;
iW=iW+1;
%See difference between initial guessed influence and local response space
%gradient
E1=Infls-Infl_f;
%Adjust Influence guess by error
Infl=Infl-1/20*E1;
%Save cumulative error
Esave(iW)=sum(E1.^2);
%Backpropagate error to determine new standard deviation guess.
StdGuess(iW)=abs(StdGuess(iW-1)-.1*(StdGuess(iW-1)-(DiffOutputs(iW)-Mean)));
end
%Save Results
True_STD1=Std_Record(iW);
w_true_z=Infl./StdGuess(iW);
Infl_z=Infl;
Efinal(iCC)=Esave(iW);
%Adjust new guessed influence by stochastic gradient descent of error
%between the influence that produces weights with the correct local gradient
%and the guessed influence
Infl=Infl_f-(Infl_f-Infl_z).*(1-Rc.^2)
%Determine Standard Deviation For fixed, new influence guess
iW=1;
%Clear as necessary
clear w
clear A DiffOutputs StdGuess
%Initialize with new guesses.
DiffOutputs(1)=gC;
StdGuess(1)=Mean;
Std_Record(1)=StdGuess(1);
%While system is above this threshold…
while abs(DiffOutputs(iW)-Mean-Std_Record(iW))>Thresh
%Based on guessed influence and guessed STD, generate vector of guessed
%input weights
w=Infl./StdGuess(iW);
%Apply guessed input weight vector to inputs
In=InputMean+Inputstd.*(w(:));
%Run Model
GcFun
%Save Outputs
Count=Count+1;
DiffOutputs(iW+1)=gC;
%Save Outputs
Std_Record(iW+1)=StdGuess(iW);
w_record(:,iW+1)=w;
iW=iW+1;
%Save Error of this guess
Esave(iW)=sum(E1.^2);
%Backpropagate Error to determine new standard deviation guess
StdGuess(iW)=abs(StdGuess(iW-1)-.1*(StdGuess(iW-1)-(DiffOutputs(iW)-Mean)));
end
%We’ve converged on the correct standard deviation
True_STD(iCC)=StdGuess(iW);
end
Appendix F: Temporal Sustainability Index MATLAB Code
%Temporal Sustainability Metrics
clc
clear all
%Import Results of Uncertainty Propagation Algorithm for Microbial
%Contamination
stdarray=xlsread(‘Array.xlsx’,’Sheet1′);
meanarray=xlsread(‘Array.xlsx’,’Sheet2′);
%Assign Safety Limit
LimitB=6.2e4;
iL=1
%Create Vector of various different Limits
LimitVector=(LimitB-LimitB/2:LimitB/15:LimitB+LimitB/2)/283.168;
%Evaluate metrics for different Limits
for Limit=LimitB-LimitB/2:LimitB/15:LimitB+LimitB/2;
%Evaluate Temporal Evolution of Risk
for it=1:length(stdarray)
%Extract the mean and standard deviation of each component at this
%point in time
mean_at_time=meanarray(it,:);
std_at_time=stdarray(it,:);
%Combine mean and standard deviation for all system components
m=sum(mean_at_time);
s=sqrt(sum(std_at_time.^2));
%For stabilitiy.
if s<0
Prob(it)=0;
else
%Determine parameters of lognormal distribution
mu = log((m^2)/sqrt(s^2+m^2));
sigma = sqrt(log(s^2/(m^2)+1));
%Determine probability of failure through cumulative lognormal
%distribution.
Prob(it)=1-logncdf(Limit,mu,sigma);
end
%Save Results
mean_sum(it)=m;
std_sum(it)=s;
end
iR=1;
%Determine for different Resilience Limits of acceptable risk
for Res_Limit=.01:.05:1;
%Vectorize Resilience Limit
ResLine(1:length(Prob))=Res_Limit;
%Find Intersect of temporal evolution of risk curve with acceptable
%risk line
IntersectX=InterX([1:length(Prob);ResLine],[1:length(Prob);Prob]);
%For Stability, make zero if empty
if isempty(IntersectX)==1
Resilience(iL,iR)=0;
Integral(iL,iR)=0;
else
%Clear as necessary
clear durations
inX=1;
%Determine duration of unacceptable probability of failure
for interii=2:2:length(IntersectX)
%An exceedance will need to have two intersects
durations(inX)=IntersectX(1,interii)-IntersectX(1,interii-1);
inX=inX+1;
end
%Determine Intersection area above curve
IntersectionsRes(iL,iR,:)=[int32(IntersectX(1,1)),int32(IntersectX(1,length(IntersectX)))];
%Determine resilience for this acceptable risk and this contamination limit
Resilience(iL,iR)=sum(durations);
%Check the Integral
Integral(iL,iR)=trapz(Prob(int32(IntersectX(1,1)):int32(IntersectX(1,length(IntersectX))))./100);
end
iR=iR+1;
end
%Vulnerability is the maximum probability of failure in respect to a
%certain threshold
Vulnerability(iL)=max(Prob);
%Save Probability
ProbMatrix(:,iL)=Prob;
iL=iL+1
end
iL=2;
iR=2;
%Convert Resilience and Vulnerability into format that can be used to
%determine sustainability index
for it=1:length(ProbMatrix(:,iL))
if it<=IntersectionsRes(iL,iR,1)
Res_time(it)=Resilience(iL,iR);
elseif it>IntersectionsRes(iL,iR,1) && it<=IntersectionsRes(iL,iR,2)
Res_time(it)=Resilience(iL,iR)-(it-IntersectionsRes(iL,iR,1));
else
Res_time(it)=0;
end
Vul_time(it)=max(ProbMatrix(1:it,iL));
end
%Determine sustainability index
SI=(1-ProbMatrix(:,iL)).*(max(Res_time)-Res_time(:))/100.*(1-Vul_time(:));
%Adjust sustainability index
SI_true=transpose(1./(.01.*Resilience(:,3))).*(1-Vulnerability)*(1-Prob(10));
hold on
%Plot results in CFU/100mL
plot((1:length(stdarray))/100,stdarray/283.168)
plot((1:length(std_sum))/100,std_sum/283.168,’k’)
xlabel(‘Time (Hour)’)
ylabel(‘Concentration (CFU/100mL)’)
title(‘Concentration Standard Deviation Curve (7/11/06)’)
figure
hold on
plot((1:length(meanarray))/100,meanarray/283.168)
plot((1:length(mean_sum))/100,mean_sum/283.168,’k’)
xlabel(‘Time (Hour)’)
ylabel(‘Concentration (CFU/100mL)’)
title(‘Concentration Mean Curve (7/11/06)’)
data=[0 0 4500 13500 60000 50000 20000 0];
time=[’12:00′ ’15:00′ ’16:00′ ’19:00′ ‘1:00’ ‘4:00’ ‘7:00′ ’10:00’];
hour=[1 4 5 8 14 17 20 23];
data_time=hour/.01;
figure
plot(.01:.05:1,transpose(1./(.01.*Resilience(1:3,:))))
