I have been reading Mixed Effects Models and Extensions in Ecology with R by Zuur et al. where departures from iid errors (heteroscedacity and/or correlation) in linear regression (and glm) are dealt with using
gls() in the
nlme package and specifically the process advocated is to set up an error correlation structure and/or weights for the model explicitly based on 1) findings from examining a residual plot or ACF, 2) conducting a likelihood ratio test using nested models or 3) AIC for non-nested models. In this manner, the analyst is able to iterate on the model selection process until the (normalized) residuals appear to satisfy the assumptions of the model.
This approach contrasts with using a robust standard error such as this blog post http://rforpublichealth.blogspot.com/2014/10/easy-clustered-standard-errors-in-r.html
My question is if both approaches are valid and simply represent different approaches or if one of the approaches is considered superior?