David Winsemius <dwinsemius <at> comcast.net> writes:
> > It's really not that difficult to get the variance covariance matrix. > What is not so clear is why you think differential weighting of a set > that has a perfect fit should give meaningfully different results than > a fit that has no weights. Again, David, what I have in mind is: since there are errors or uncertainties in the response variables (despite the perfect collinearity of the data), which I assume are Gaussian, if I make a large enough number of simulations of four response values, there will undoubtedly be a dispersion in the best fit intercept and slope obtained from a usual unweighted least squares procedure, right? Then, if I calculate the arithmetic mean of these simulated intercept and slope, I would certainly check that they would be 0 and 2, respectively. However, and THAT IS THE POINT, there will also be a standard deviation associated with each one of these two coefficients, right??, and that is what I would assign as the measure of uncertainty in the estimation of the coefficients. This is not, as Dalgaard has called attention to, what the simple command summary(lm(y~x,weights=1/err^2)) provides in its Std. Error. However, as Dalgaard also recalled, the command summary(glm(y~x,family=gaussian,weights=1/err^2),dispersion=1) does provide Std. Errors in the coefficients which look plausible (at least to me) and, at any rate, which do coincide with results from other packages (Numerical Recipes, ROOT and possibly GSL...) > > ?lm > ?vcov > > > y <- c(2,4,6,8) # response vect > > fit_mod <- lm(y~x,weights=1/error^2) > Error in eval(expr, envir, enclos) : object 'error' not found > > error <- c(0.3,0.3,0.3,0.3) > > fit_mod <- lm(y~x,weights=1/error^2) > > vcov(fit_mod) > (Intercept) x > (Intercept) 2.396165e-30 -7.987217e-31 > x -7.987217e-31 3.194887e-31 > > Numerically those are effectively zero. > > > fit_mod <- lm(y~x) > > vcov(fit_mod) > (Intercept) x > (Intercept) 0 0 > x 0 0 > ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
