Date: April 15, 2013 Memo 6-7.1

To: Cuc Schroeder, EPA

From: John Rogers

Subject: Serial correlations for Steam Electric with and without adjustment for baseline values

[EP-C-10-023 WA 2-6 (Steam Electric)] 

Revisions: This memo reflects wording changes and clarifications describing the correlation 

estimates originally delivered to EPA in Memo 6-4 dated September 27, 2012. The most recent 

change revised the tables to show only the correlations used in the development the final effluent 

limitations. 

1. Introduction

Under EPA Contract EP-C-10-023 WA 2-6, EPA has asked Westat to calculate the 1-day lag 

serial correlations to be used when calculating limits for steam electric plants and to calculate the 

limits for various technology options defined by EPA. This memo describes the process for 

calculating the correlations and provides estimates of the serial correlation using 1) the data, 

adjusted by a set of baseline values specified by EPA that were used for calculating the proposed 

limits for the Steam Electric Effluent Guidelines rulemaking; and 2) all valid reported results, i.e., 

the data collected by EPA or provided by the plants with possible non-detects but without 

baseline substitution. The autocorrelations obtained under these two scenarios were used to 

calculate two sets of limits (with and without adjustment for baseline values).

EPA used the baseline values in the statistical analysis to account for the possibility that certain 

detected results may be at a lower concentration than can be reliably achieved by well-operated 

laboratories. For calculating the limits from the data, both detected and non-detected values that 

are below the baseline are set equal to the baseline and treated as not detected. In this memo, the 

term “non-detect” or “not detected” refers to 1) values below the detection limit in the original 

reported results or 2) values below the detection limit in the original reported results or below the 

selected baseline, depending on the context. 

Table 1 shows the baseline values established by EPA for the calculation of the serial correlations 

and limits. 
-2-

Table 1 Baseline values by pollutant

Pollutant Baseline value Unit

Arsenic 2 ug/L

Selenium 5 ug/L

Mercury 0.5 ng/L

Nitrate Nitrite as N 0.05 mg/L

Total Dissolved Solids 10 mg/L

For some pollutants and treatment technology options for FGD and gasification wastewater, the 

autocorrelations could not be calculated due to a limited amount of data. As a result, serial 

correlations could be calculated for only arsenic, selenium, and mercury for the biological 

treatment technology option for FGD wastewater. For nitrate-nitrite as N for this technology 

option, since there were only 8 observations available, autocorrelation could not be calculated. For 

the biological technology option for FGD wastewater, EPA is proposing limits for arsenic, 

mercury, nitrate-nitrite as N, and selenium. For arsenic and mercury, EPA transferred the limits 

from the chemical precipitation treatment system to the biological treatment for FGD wastewater. 

As a result, this memo only presents the serial correlation estimates for selenium for the biological 

treatment option. Of the baseline values considered by EPA, only those for arsenic and selenium 

affect the calculation of serial correlations and limits. The baseline for mercury and other 

pollutants are less than the smallest observed data value. This memo presents the serial 

correlations for selenium adjusted for the baseline values shown in Table 1 and using all valid 

reported results (without baseline adjustment).

2. Data Processing

The data processing and analysis was performed using SAS 9.2. Westat received the data from 

EPA in the SE_Sampling_Data_10_1_12.xls Excel file. The data were imported into SAS and 

processed by: 1) dropping values to be excluded from the calculations; 2) adjusting the 

concentrations for the baseline values; 3) creating some categorical variables to facilitate the data 

processing, 4) averaging duplicate samples collected on the same day by EPA, 5) averaging the 

results of two samples from the same day, one collected for EPA and one collected by the plant; 

and 6) calculating the log-transformed concentration.

3. Calculation of 1-day-lag serial correlation

For a series of daily effluent measurements with all detected values, the calculation of the 1-day lag 

serial correlation assumes the log-transformed data can be described by an autoregressive model 

with one correlation term for 1-day lag, an AR(1) model. The model assumes the log-transformed 

data have a normal distribution. The AR(1) model has been used by EPA in the past to describe 

the serial correlation in effluent series. A preliminary analysis suggests that the AR(1) model is also 

a reasonable choice for the steam electric data. 

In general, the serial correlation term can be estimated using time series ARIMA models. 

However, for the steam electric data, the estimation is complicated by small sample sizes, gaps in 
-3-

the series, and non-detected values. With small sample sizes (particularly less than 30) the 

algorithms can have difficulty converging, the estimates will be particularly imprecise, and the 

estimates will, on average, be closer to zero than the true correlation. The presence of non-detects 

in the series is particularly problematic. In many analyses, the non-detects are replaced by a 

substitute value, such as the detection limit (or baseline) or the detection limit over two. The 

presence of non-detects and the choice of the substitute value can bias the correlation estimates. 

EPA generally substitutes the detection limit (or baseline) for the non-detected values. If there are 

several different detection limits, the estimated serial correlation may be more closely related to 

the serial correlation of the detection limits than to the correlation of the concentrations being 

measured.

The serial correlation is usually obtained by finding the parameter estimate that maximizes the 

likelihood of the data. The standard ARIMA models assume a likelihood function that is correct 

for serially correlated data with only detected values. Alternatively, the data can be analyzed using 

survival analysis using the SAS LIFEREG procedure. When using survival analysis, different 

likelihood functions are used for the detects and the non-detects where the likelihood function for 

the non-detects depends on the detection limit. Survival analysis estimates the parameters of the 

hypothetical distribution that would have been observed if the field samples had been analyzed 

using a more precise analytical method such that no observations were below the detection limit. 

Although survival analysis uses a likelihood function that is correct for data with non-detects, it 

does not account for or estimate the serial correlation. Westat modified an algorithm developed in 

other work which combines the likelihood function for serially correlated data with the likelihood 

function for non-detects. Based on simulations, the Westat procedure provides better estimates of 

the serial correlation for the concentration series than when using the ARIMA model with 

substitute values for the non-detects. The Westat procedure was implemented using a Bayesian 

approach as described below.

Note that the correlation estimated by the Westat algorithm is the serial correlation for the series 

that we would have observed if the samples had been analyzed using a laboratory technique with 

low detection limit such that all concentrations were detected (call this the underlying series). The 

adjustment for non-detects is achieved by using a different likelihood function for the detected 

and the non-detected values. Call the resulting estimated correlation value for the underlying series

Rho. However, this is not the serial correlation that is needed for the calculating the effluent 

limits. The limit calculations use the serial correlation of the series with values below the baseline 

replaced by the baseline and non-detects above the baseline replaced by the detection limit (call 

this the baseline adjusted data). Call the serial correlation for the baseline adjusted data the 

baseline adjusted Rho. To estimate the baseline adjusted Rho we can either: 1) estimate Rho using 

the Westat macro applied to all valid reported results and adjust the estimate to account for the 

presence of non-detects (as described below in Section 6); or 2) calculate the baseline adjusted 

Rho from the baseline adjusted data. Although the second approach is relatively straight forward, 

the process of setting values below the baseline equal to the baseline reduces the information in 

the data and as a result the correlation estimates using the second approach are less precise than 

from the first approach (this was verified using simulations). In addition, the estimates can be 

affected by possible serial correlation in the detection limits. The first approach avoids these 

problems. Although the correlations calculated both ways are similar, the final correlations are 

based on estimates of Rho which are subsequently adjusted for the presence of non-detects, as 

described below. 
-4-

4. The Bayesian approach for estimating Rho 

The more common and traditional statistical analysis procedures, called frequentist procedures, 

assume the parameters are fixed and the data are subject to random error. The resulting parameter 

provides an estimate of the unknown parameter and the 95% confidence interval provides a range 

which is expected to cover the true parameter value 95% of the time. In contrast, Bayesian 

procedures assume the data are fixed and the parameter estimates have a distribution representing 

the possible values of the parameter given the data and any prior information about the parameter 

(the Bayesian model we used assumes essentially no prior information, using a very diffuse prior 

distribution). For the non-detects, the Bayesian model uses a substitute value to aid in estimating 

the serial correlation but not for the likelihood for the non-detects. Calculation of the parameter 

distributions for the Bayesian model involves integrating the likelihood function, which can be 

very difficult. As a result, the distributions are generally approximated using simulation. The mean 

and standard deviation of the simulated values are often used to summarize our knowledge of the 

parameter value. In practice, the frequentist and Bayesian approaches provide similar results. For 

this problem, the Bayesian approach was easier to adapt than a frequentist approach for obtaining 

estimates of the serial correlation when some values are below the detection limit. Westat created 

a macro to estimate the serial correlation in data with non-detects and with gaps in the data which 

used Bayesian methods implemented using the SAS MCMC procedure. To minimize the 

estimation error, the correlations are based on 200,000 simulated values. 

5. Simulations to evaluate the Bayesian approach for estimating Rho 

The validity of this approach was assessed by simulating data with different characteristics and 

comparing two analysis approaches: 

1) Analyze the data using the SAS ARIMA procedure with a substitute value for the nondetects (or baseline); or

2) Analyze the data using the Westat Bayesian approach to estimate Rho and adjusting Rho 

for the presence of values below the detection limit.

The simulated series of daily data used the following parameters:

• Set the substitute value to the median of the expected distribution of the non-detects;

• 0%, 50%, or 80% of the daily values were randomly missing

• 25, 50, or 200 non-missing observations;

• Rho = 0, .3, .6, or .9; and

• On average, 0%, 25%, 50%, or 75% of the non-missing observations were below the 

detection limit.

The results are based on multiple simulated datasets for each combination of parameters above 

(10 datasets when simulating 200 non-missing observations, 40 when simulating 50 non-missing 

observations, and 80 simulated datasets when simulating 25 non-missing observations). These two 

approaches were compared using the bias in the estimated Rho (estimated Rho minus the true 

Rho used to generate the data) and the bias in the estimated standard error. If the standard error 

estimate from the ARIMA or Bayesian procedure is a good estimate of the true uncertainty in the 

estimate, the ratio of the pooled standard error across simulated datasets over the standard 
-5-

deviation of the estimates from the simulated datasets would be 1.0. The log of this ratio was used 

as an estimate of the bias in the estimated standard error.

In general, for both approaches:

1) The procedures performed well with little bias with 200 non-missing observations and no 

or few non-detects. 

2) The estimated standard error was more likely to underestimate the true error with high 

correlations and, for the ARIMA model, small sample sizes. 

3) The estimated Rho was more likely to underestimate the true Rho for smaller sample 

sizes, higher proportions of non-detects, and higher correlations. 

Across all simulated conditions, the Bayesian approach provided a similar or less biased estimate 

of Rho and a more appropriate (or even conservative) estimate of the standard error. Across all 

simulations, the average bias provides a good comparison of the approaches. The average bias 

across all simulations was -0.036 for the Bayesian approach and -0.084 for the ARIMA approach. 

The average bias of the standard error was 0.027 for the Bayesian approach (i.e., standard errors 

were generally higher (conservative) by 2.7%) and -0.184 for the ARIMA approach (standard 

errors were generally low by 18%).

The estimates from the Bayesian approach are essentially unbiased for data similar to the available 

data, i.e., 50 non-missing values, of which half were below the detection limit, spread across 100 

days, with true Rho of 0.6. In contrast, using the ARIMA model, the estimates of Rho are low by 

0.09 and the standard errors are low by roughly 28%.

The objective of the simulations was to verify that the Bayesian approach provides estimates of 

Rho that are comparable to or better than a simple ARIMA model for estimating Rho. Based on 

the simulations, we expect the Bayesian estimates to be more precise and less biased. Although the 

results depend somewhat on the substitute value used for non-detects, we expect the Bayesian 

approach to perform at least as well as the ARIMA approach. We also expect the Bayesian 

estimates to be less biased because the Bayesian procedure explicitly adjusts for the presence of 

non-detects. As a result, the Bayesian approach was used to estimate the serial correlations used 

when calculating the limits.

6. Calculation of the baseline adjusted Rho

The baseline adjusted Rho was calculated from the estimate of Rho and the proportion of nondetects (P) by:

1) Estimate Rho by applying the Bayesian algorithm to all valid reported results to estimate 

the serial correlation for the underlying series.

2) Calculating the proportion of non-detects for a selected baseline (P) from the baseline 

adjusted data;

3) Simulating one million sequential daily observations with a standard normal distribution 

and serial correlation Rho;

4) Setting values below the Pth percentile to the Pth percentile of the distribution; and

5) Setting the baseline adjusted Rho equal to the 1-day lag correlation of the resulting series, 

calculated using the SAS ARIMA procedure.
-6-

7. Correlations Used for the Limits

Separate correlations were calculated for each combination of treatment technology option, 

pollutant, plant, and baseline for FGD and gasification wastewaters. However, for many of these 

combinations, there were very few observations (e.g., 8) or very few observations above the 

detection limit. Estimates based on such data would be very imprecise, or in some cases, reliable 

estimates could not be obtained due to convergence problems. As a result, Rho was estimated 

only for datasets with a reasonable number of detected values and total number of observations. 

For the estimates that were calculated, there were at least 20 detected values and at least 38 values 

overall. There was enough data to calculate the serial correlation for only biological treatment for 

FGD wastewater.

For biological treatment for FGD wastewater, Table 2 shows the estimated serial correlations and 

standard errors (in parentheses). In Table 2, Rho is the estimated correlation for the underlying 

series (the series that we would have observed if the samples had been analyzed using a laboratory 

technique with a low detection limit such that all concentrations were detected).

Table 2. Correlation estimates for the biological treatment for FGD wastewater

Analyte Plant Name N N 

Detected Rho (±SE)

Selenium Allen 63 20 0.8552 (±0.0501)

Belews Creek 77 36 -0.1902 (±0.2946)

Table 3 shows the estimated serial correlation for the baseline adjusted data. In Table 3, Baseline

is the assumed baseline value, N detected is the number of detected values after adjusting for the 

baseline. Baseline Adjusted Rho is the estimated correlation for the baseline adjusted data. For a 

treatment option and pollutant, if Baseline Adjusted Rho is estimated for multiple plants, EPA 

uses the mean correlation across the plants on the assumption that the correlation can vary 

between plants and the mean is a better estimate of the correlation that might be found at plants 

covered by the regulations. The mean is shown in the last column (Mean Baseline Adjusted 

Rho) and is the correlation used for the limits. A baseline equal to zero is equivalent to not using a 

baseline. 

Table 3. Baseline adjusted correlation estimates for the biological treatment technology option for 

FGD wastewater

Pollutant Plant Name N Baseline

N detected 

(After 

Baseline 

Adjustment)

Baseline 

Adjusted

Rho

Mean

Baseline 

Adjusted

Rho

Selenium

ug/L

Allen 63 0 20 0.7835 0.329 Belews Creek 77 36 -0.1254

Allen 63 54 0.6815 0.291 Belews Creek 77 26 -0.0997