[R] Uncertainty propagation

[Ben Bolker]

Maayt <m.lupker <at> hotmail.com> writes:

> I linearized my power relations en fitted them with a linear rlm() function.
> When I re-sample my pairs from a bivariate normal distribution for my power
> law what transformation do I need to apply a transformation to my covariance
> (vcov) matrix to get back from my linearized regression to my power law
> "space".
> 
> Thanks

   rlm() does return a fitted model object that 'inherits from' (is
a variant/superset of) the 'lm' class, therefore vcov(modelfit) should
work (but see caution below).  You should 'unpack' the results in the opposite
direction from your modeling -- simulate on the linearized scale, then
invert the transformation you used in order to get the curves back to your 
'power law space'.
   The caution is that without digging into the details of rlm() [and
reading the appropriate section of Venables and Ripley], I don't know 
whether the vcov() matrix based on robust regression will preserve the
non-Gaussian characteristics of your data ... you may find when you do
the simulations that they do *not* capture the variance of your data
appropriately, because the robust part of 'robust linear modeling' 
deliberately downweights the effects of outliers.
   You may find that your results look plausible anyway.
   If not, this begins to turn (for me anyway) into a non-trivial
problem.  One possibility (although more time-consuming) would be
to (nonparametrically) bootstrap the data, and generate a predicted
curve for each bootstrap sample -- then use the envelope of these
bootstrapped curves to characterize the uncertainty (in general this
would be a little bit more robust/general/parsimonious
than the 'prediction interval' approach I've suggested, although more
computationally intensive and slightly harder to set up).