Mark Delorey F. Jay Breidt Colorado State University

1. This research is funded by U.S.EPA – Science To Achieve Results (STAR) Program Cooperative Agreement # CR - 829095 Distribution Function Estimation in Small Areas for Aquatic ResourcesSpatial Ensemble Estimates of Temporal Trendsin Acid Neutralizing Capacity Mark Delorey F. Jay Breidt Colorado State University

2. Project Funding • The work reported here was developed under the STAR Research Assistance Agreement CR-829095 awarded by the U.S. Environmental Protection Agency (EPA) to Colorado State University. This presentation has not been formally reviewed by EPA.  The views expressed here are solely those of the presenter and STARMAP, the Program he represents. EPA does not endorse any products or commercial services mentioned in this presentation.

3. Outline • Statement of the problem: How to get a set of estimates that are good for multiple inferences of acid trends in watersheds? • Hierarchical model and Bayesian inference • Constrained Bayes estimators • adjusting the variance of the estimators • Conditional auto-regressive (CAR) model • introducing spatial correlation • Constrained Bayes with CAR • Summary

4. The Problem • Examine acid neutralizing capacity • Supply of acids from atmospheric deposition andwatershed processes exceeds buffering capacity • Surface waters are acidic if ANC < 0 • Temporal trends in ANC within watersheds (8-digit HUC’s) • characterize the spatial ensemble of trends • make a map, construct a histogram, plot an empirical distribution function

5. Data Set • 86 HUC’s in Mid-Atlantic Highlands • ANC in at least two years from 1993–1998 • HUC-level covariates: • area • average elevation • average slope, max slope • percents agriculture, urban, and forest • spatial coordinates

6. Region of Study

7. Locations of Sites

8. Small Area Estimation • Probability sample across region • regional-level inferences are model-free • samples are not sufficiently dense in small watersheds (HUC-8) • need to incorporate auxiliary information through model • Two standard types of small area models (Rao, 2003) • area-level: watersheds • unit-level: site within watershed

9. Two Inferential Goals • Interested in estimating individual HUC-specific slopes • Also interested in ensemble:spatially-indexed true values:spatially-indexed estimates: • subgroup analysis: what proportion of HUC’s have ANC increasing over time? • “empirical” distribution function (edf):

10. Deconvolution Approach • Treat this as measurement error problem: • Deconvolve: • parametric: assume F in parametric class • semi-parametric: assume F well-approximated within class (like splines, normal mixtures) • non-parametric: assume EF[ei] is smooth • Not so appropriate for heteroskedastic measurements, explanatory variables, two inferential goals

11. Hierarchical Area-Level Model • Extend model specification by describing parameter uncertainty: • Prior specification:

12. Bayesian Inference • Individual estimates: use posterior meanswhere • Do Bayes estimates yield a good ensemble estimate? • use edf of Bayes estimates to estimate F? • No: Bayes estimates are “over-shrunk” • too little variability to give good representation of edf (Louis 1984, Ghosh 1992)

13. Adjusted Shrinkage • Posterior means not good for both individual and ensemble estimates • Improve by reducing shrinkage • sample mean of Bayes estimates already matches posterior mean of • adjust shrinkage so that sample variance of estimates matches posterior variance of true values • Louis (1984), Ghosh (1992) • Cressie and Stern (1991)

14. Constrained Bayes Estimates • Compute the scalars • Form the constrained Bayes (CB) estimates aswhere

15. Shrinkage Comparisons for the Slope Ensemble

16. Numerical Illustration • Compare edf’s of estimates to posterior mean of F: • Comparison of ensemble estimates at selected quantiles:

17. Estimated EDF’s of the Slope Ensemble CB Posterior Mean Bayes

18. Conditional Auto Regressive (CAR) Model • Letwhere  is an unknown coefficient vector, C = (cij) represents the adjacency matrix,  is a parameter measuring spatial dependence,  is a known diagonal matrix of scaling factors for the variance in each HUC, and  is an unknown parameter. • Adjacency matrix C can reflect watershed structure

19. HUC Structure • First level (2-digit) divides U.S. into 21 major geographic regions • Second level (4-digit) identifies area drained by a river system, closed basin, or coastal drainage area • Third level (6-digit) creates accounting units of surface drainage basins or combination of basins • Fourth level (8-digit) distinguishes parts of drainage basins and unique hydrologic features

20. Neighborhood Structure • All watersheds within the same HUC-6 region were considered part of same neighborhood • No spatial relationship among HUC-4 regions or HUC-2 regions considered at this point

21. Model Specifications • Adjacency matrix: •  = Im •  fixed at three different levels: 0.09, 0.01, 0

22. Constrained Bayes with CAR • Again, compute H1 and H2:where • Then,

23. Comparison of Estimated EDF’s

24. Spatial Structure

25. Summary • In Bayesian context, posterior means are overshrunk; in order to obtain estimates appropriate for ensemble, need to adjust • In CAR model, value of  does have an effect on edf, but not large effect, possibly due to CB calculation • Edf of CB estimates in CAR shifts more mass towards positive values • Contour plot indicates that trend slopes of ANC are smoothed within HUCs

26. Further Work • Restrict to acid-sensitive waters • Combine probability and convenience samples • Other covariates? • deposition maps/trends from CASTNet? • Modify spatial structure • Site-level model? • useful sub-watershed covariates? • spatial scales: HUC to HUC, site to site • more concern with design, normality assumptions