G'day Daniel, >>>>> "DR" == daniel <[EMAIL PROTECTED]> writes:
DR> First: I think coefficients from lm(Employed~.,data=longley) DR> should be equal coefficients from DR> lm.ridge(Employed~.,data=longley, lambda=0) why it does not DR> happen? Which version of R and which version of MASS are you using? > lm(Employed~.,data=longley) Call: lm(formula = Employed ~ ., data = longley) Coefficients: (Intercept) GNP.deflator GNP Unemployed Armed.Forces -3.482e+03 1.506e-02 -3.582e-02 -2.020e-02 -1.033e-02 Population Year -5.110e-02 1.829e+00 > lm.ridge(Employed~.,data=longley, lambda=0) GNP.deflator GNP Unemployed Armed.Forces -3.482259e+03 1.506187e-02 -3.581918e-02 -2.020230e-02 -1.033227e-02 Population Year -5.110411e-02 1.829151e+00 These coefficients look pretty identical to me, except that they are printed to different numbers of significant digits. In fact, the following shows that they are identical (upto numerical precision): > fm1 <- lm(Employed~.,data=longley) > fm2 <- lm.ridge(Employed~.,data=longley, lambda=0) > coef2 <- print(fm2) GNP.deflator GNP Unemployed Armed.Forces -3.482259e+03 1.506187e-02 -3.581918e-02 -2.020230e-02 -1.033227e-02 Population Year -5.110411e-02 1.829151e+00 > max(abs(coef(fm1)-coef2)) [1] 7.275958e-12 DR> Second: if I have for example Ridge<-lm.ridge(Employed~., DR> data=longley, lambda = seq(0,0.1,0.001)), I suppose intercept DR> coefficient is defined implicit, Yes. DR> why it does not appear in Ridge$coef? If you look at the code of lm.ridge, you will see that, if an intercept is included in the model, all non-constant regressors are centered (i.e. made orthogonal to the intercept term) and scaled to have the same variance. Further more, the intercept term is typically *not* penalised. The components in Ridge$coef are the coefficients on this transformed scale. No need of including the intercept here, since it is the same for all values of lambda. If you print the model, then the ridge coefficients on the original scale are calculated, see: > getAnywhere("print.ridgelm") A single object matching 'print.ridgelm' was found It was found in the following places registered S3 method for print from namespace MASS namespace:MASS with value function (x, ...) { scaledcoef <- t(as.matrix(x$coef/x$scales)) if (x$Inter) { inter <- x$ym - scaledcoef %*% x$xm scaledcoef <- cbind(Intercept = inter, scaledcoef) } print(drop(scaledcoef), ...) } <environment: namespace:MASS> DR> Third: I suppose that if I define DR> 1) y<-longley$Employed DR> 2) X<-as.matrix(cbind(1,Longley[,1:6]) DR> 3) I = identity matrix the DR> following should be true: Coef=(X'X+kI)^(-1) X'y No, as noted above, the intercept term is usually not penalised. DR> and if a take k=Ridge$kHKV, Coef should be approx equal to DR> Ridge$Coef[near value of kHKV] No, as noted above the estimates in the "coef" component of an object returned by lm.ridge are the coefficients on a different scale. DR> and it does not seem to happen, why? Because the intercept is not penalised by lm.ridge and the non-constant columns of the design matrix are rescaled; hence the returned coefficients are on another scale. DR> Any help, suggestion or orientation? HTH. Cheers, Berwin ______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html