R help - Aug 2009 - ridge regression

Dear all,

I considered an ordinary ridge regression problem. I followed three
different ways:
1. estimate beta without any standardization
2. estimate standardized beta (standardizing X and y) and then again convert
back  
3. estimate beta using lm.ridge() function


X<-matrix(c(1,2,9,3,2,4,7,2,3,5,9,1),4,3)
y<-t(as.matrix(cbind(2,3,4,5)))

n<-nrow(X)
p<-ncol(X)

#Without standardization
intercept <- rep(1,n)
Xn <- cbind(intercept, X)
K<-diag(c(0,rep(1.5,p)))
beta1 <- solve(t(Xn)%*%Xn+K)%*%t(Xn)%*%y
beta1

#with standardization
ys<-scale(y)
Xs<-scale(X)
K<-diag(1.5,p)
bs <- solve(t(Xs)%*%Xs+K)%*%t(Xs)%*%ys 
b<- sd(y)*(bs/sd(X))
intercept <- mean(y)-sum(as.matrix(colMeans(X))*b)
beta2<-rbind(intercept,b)
beta2

#Using lm.ridge function of MASS package
beta3<-lm.ridge(y~X,lambda=1.5)

I'm getting three different outputs for using the above three different
approaches, but which would been the same for all. 
> beta1
                [,1]
intercept  3.4007944
           0.3977462
           0.2082025
          -0.4829115
> beta2
                [,1]
intercept  3.3399855
           0.1639469
           0.0262021
          -0.1228987
> beta3
                     X1          X2          X3 
 3.35158977  0.19460958  0.03152778 -0.15546775 

It will be very helpful to me if anybody can help me regarding why the
outputs are coming different. 

regards.
-- 
View this message in context:
http://www.nabble.com/ridge-regression-tp25045020p25045020.html
Sent from the R help mailing list archive at Nabble.com.

If you didn't post anonymously I would have made a suggestion.  Full 
names and affiliations should be given.

Frank


spime wrote:
> Dear all,
> 
> I considered an ordinary ridge regression problem. I followed three
> different ways:
> 1. estimate beta without any standardization
> 2. estimate standardized beta (standardizing X and y) and then again
convert
> back  
> 3. estimate beta using lm.ridge() function
> 
> 
> X<-matrix(c(1,2,9,3,2,4,7,2,3,5,9,1),4,3)
> y<-t(as.matrix(cbind(2,3,4,5)))
> 
> n<-nrow(X)
> p<-ncol(X)
> 
> #Without standardization
> intercept <- rep(1,n)
> Xn <- cbind(intercept, X)
> K<-diag(c(0,rep(1.5,p)))
> beta1 <- solve(t(Xn)%*%Xn+K)%*%t(Xn)%*%y
> beta1
> 
> #with standardization
> ys<-scale(y)
> Xs<-scale(X)
> K<-diag(1.5,p)
> bs <- solve(t(Xs)%*%Xs+K)%*%t(Xs)%*%ys 
> b<- sd(y)*(bs/sd(X))
> intercept <- mean(y)-sum(as.matrix(colMeans(X))*b)
> beta2<-rbind(intercept,b)
> beta2
> 
> #Using lm.ridge function of MASS package
> beta3<-lm.ridge(y~X,lambda=1.5)
> 
> I'm getting three different outputs for using the above three different
> approaches, but which would been the same for all. 
> 
>> beta1
>                 [,1]
> intercept  3.4007944
>            0.3977462
>            0.2082025
>           -0.4829115
>> beta2
>                 [,1]
> intercept  3.3399855
>            0.1639469
>            0.0262021
>           -0.1228987
> 
>> beta3
>                      X1          X2          X3 
>  3.35158977  0.19460958  0.03152778 -0.15546775 
> 
> It will be very helpful to me if anybody can help me regarding why the
> outputs are coming different. 
> 
> regards.

-- 
Frank E Harrell Jr   Professor and Chair           School of Medicine
                      Department of Biostatistics   Vanderbilt University