320 likes | 511 Views
Bayesian and Classical Analysis of Multi-Stratum Response Surface Designs. Steven Gilmour Queen Mary, University of London Peter Goos Universiteit Antwerpen. Outline. Split-plot and other multi-stratum designs “State-of-the-art” analysis of data REML/generalized least squares Problems
E N D
Bayesian and Classical Analysis of Multi-Stratum Response Surface Designs Steven Gilmour Queen Mary, University of London Peter Goos Universiteit Antwerpen
Outline • Split-plot and other multi-stratum designs • “State-of-the-art” analysis of data • REML/generalized least squares • Problems • Estimation of variance components • Degrees of freedom • Three possible solutions • Fix value(s) of variance-component(s) • Use randomization-based estimation • Bayesian analysis
Multi-stratum designs • Randomization of treatments to experimental units is restricted in such a way that particular sets of units must receive the same level of one or more treatment factors • Includes classical orthogonal split-plot, split-split-plot, criss-cross, etc. designs (regular factorial treatment sets) • Also includes nonorthogonal designs with similar structures (irregular factorial or response surface treatment sets) • Are (nested) block designs with at least one main effect totally confounded with block effects • Often necessary when some factors are hard to change
Multi-stratumdesigns • I refer to the runs as units and the groups of units defined by the randomization restrictions as blocks, superblocks, … • Randomization is performed by randomly relabelling …, superblocks, blocks and units • Implies random effects for …, superblocks, blocks, units (error) in derived linear model • Fixed treatment effects can be modelled using usual polynomial response surface model
Freeze-dried coffee experiment • Response: amount of retained volatile compounds in freeze-dried coffee • Treatment factors: • Pressure in drying chamber (dial-controlled) • Heating temperature, Initial solids content, Slab thickness, Freezing rate (all easy to change) • 5 runs during each of 6 days • Randomization restricted so that all runs in a day have the same pressure
Model • Generalized least squares (GLS) estimation • Variance-covariance matrix Model and analysis
Model • Generalized least squares (GLS) estimation • Variance-covariance matrix Model and analysis ? ?
Variance component estimation • REML: REsidual Maximum Likelihood • Yields the same answers as ANOVA in orthogonal designs (e.g. standard split-plots) • Applicable when designs are not orthogonal (e.g. nonorthogonal split-plots) • State of the art in many disciplines • Available in many statistical software packages
Analysis • Different implementations: • Variance components allowed to be negative or not • Various methods for obtaining effective degrees of freedom • Estimates generally consistent with each other, given different implementations • Effective degrees of freedom can be inconsistent with each other • All methods can give surprising results
Freeze-dried coffee experiment • Estimates of whole-plot error variance • SAS: 0 • GenStat: 0 • R: 0.0051 • Degrees of freedom for testing linear effect of pressure (full model) • SAS proc mixed with Kenward & Roger: 9 df • SAS proc mixed with containment method: 6 df • R lme default: 3 df
Freeze-dried coffee experiment Simplified model:
Freeze-dried coffee experiment • Data are treated as if they come from a completely randomized experiment • OLS estimates are obtained • Degrees of freedom for testing linear effect of pressure are too optimistic • Upper bound for full second-order model: 3 df • Because of nonorthogonality: less than 3 df • SAS proc mixed with Kenward & Roger: 9 df • SAS proc mixed with containment method: 6 df • R lme default: 3 df
Solution II: Randomization-based analysis • Even nonorthogonal multi-stratum designs have simple orthogonal block structures (if each block/superblock/... is the same size) [Nelder, 1965] • Ignoring treatment structure, randomization-based analysis gives minimum variance unbiased estimators of variance components (pure error) • Only assumption is that treatment and unit effects are additive
Randomization-based analysis • Proposed analysis: • Use discrete treatments defined by combinations of factor levels (ignoring treatment model) • Anova gives correct estimates of variance components with correct degrees of freedom • Use these estimates to fit treatment model using GLS • Base inferences on these estimates • “Extra sums of squares” represent lack of fit • Not clear that GLS is best, but is same as with REML
Freeze-dried coffee experiment • There are 27 treatments, 3 replicated twice • 0 residual degrees of freedom for blocks • 3 residual degrees of freedom for runs • Blocks variance component cannot be estimated, unit variance badly estimated • Full polynomial model can be fitted, but no global inference is possible • Weak inference is possible for all individual parameters except main effects of pressure • This design is too small for a frequentist analysis
Solution III: Bayesian approach • Advantages: • Takes into account uncertainty in prior beliefs • Prior beliefs can be contradicted by the data • No problems determining the appropriate degrees of freedom for hypothesis tests • WinBUGS software is free
The Bayesian approach • Requires a user-specified (joint) distribution for all model parameters (, 2, 2) • Posterior marginal distributions can be used for inference about parameters • Results: • Similar to REML/GLS if data contain enough information • Similar to prior distribution if data don’t contain enough information
The Bayesian approach • Noninformative priors for r: N(0,) • Weakly informative priors for r: • Linear and interaction effects: N(0,25) • Quadratic effects: N(0,100)
The Bayesian approach • Variance components
weakly informative highly informative not informative The Bayesian approach • Variance components
Summary of results • Prior information on has little impact • Prior information on 2 not important at all • Some results strongly depend on prior information about 2 • Hard-to-change factor coefficients • Sub-plot factor interaction coefficients that are not nearly orthogonal to whole plots • Results for other coefficients insensitive to the choice of the prior for 2
Discussion • REML/GLS analysis can be misleading as it often leads to an analysis that ignores the multi-stratum nature of the design • Likelihood methods have good asymptotic properties, i.e. large numbers of units in each stratum, so should not be expected to work in small experiments • Problem is due to a lack of information in the blocks stratum • We should honestly admit that there is no information and/or provide prior information
Discussion • Randomization based analysis should always be done (in every experiment!) as a first step • Makes very few assumptions, so is much more robust than any other analysis • Provides a “reality check” • Might make extra assumptions unnecessary • Bayesian analysis can help • Prior information is taken into account • Prior information can be overruled • Depends heavily on prior assumptions, but these are clearly and honestly expressed
References • Multi-stratum response surface designs Luzia A. Trinca and Steven G. Gilmour Technometrics, 2001 • A split-unit response surface design for improving aroma retention in freeze dried coffee Steven G. Gilmour, J. Mauricio Pardo, Luzia A. Trinca, K. Niranjan and Don Mottram Proceedings of the 6th European Conference on Food-Industry and Statistics, 2000 • Analysis of data from unbalanced multi-stratum designs Steven G. Gilmour and Peter Goos Submitted