2
$\begingroup$

I was able to reproduce table 3.1 from ESL. However, when I tried to reproduce table 3.2, my estimated coefficients were way off (shown below):

         [,1]
 [1,]  0.4292
 [2,]  0.5765
 [3,]  0.6140
 [4,] -0.0190
 [5,]  0.1448
 [6,]  0.7372
 [7,] -0.2063
 [8,] -0.0295
 [9,]  0.0095

The results from table 3.2 in ESL are as follow: enter image description here

I have attached my code here:

res <- read.table("prostate.data", sep = "")
fix(res)
XTraining = subset(res, train)
XTesting = subset(res, train == FALSE)


nrow = dim( XTraining )[1]
p = dim( XTraining )[2] - 1 # the last column is the response 

D = XTraining[,1:p-1] # get the predictor data

# This gives the raw data for Table 3.1 from the book:
# 
print(cor(D),digts=3)

library(xtable)
xtable( cor(D), caption="Duplication of the values from Table 3.1 from the book", digits=3 )

#
# Duplicate Table 3.2 from the book 
#

# Append a column of ones:
# 
Dp = cbind( matrix(1,nrow,1), as.matrix( D ) )

lpsa = XTraining[,p]

library(MASS)

betaHat = ginv( t(Dp) %*% Dp ) %*% t(Dp) %*% as.matrix(lpsa)

# this is basically the first column in Table 3.2:
#
print('first column: beta estimates')
print(betaHat,digits=2)

Does anyone have similar issues ? I would like to know what went wrong with my output.

$\endgroup$
6
$\begingroup$

The problem is you have not scaled the data as described in the information file:

http://statweb.stanford.edu/~tibs/ElemStatLearn/datasets/prostate.info.txt

The features must first be scaled to have mean zero and variance 96 (=n) before the analyses in Tables 3.1 and beyond. That is, if x is the 96 by 8 matrix of features, we compute xp <- scale(x,TRUE,TRUE)

This gives the correct answer:

X <- read.table("prostate.data", sep = "")
train <- X[,10]
y <- X[,9]
X <- X[,-(9:10)]

X.scaled <- scale(X,TRUE,TRUE)

betaHat <- lm(y[train]~X.scaled[train,])$coef
print(betaHat)

         (Intercept)  X.scaled[train, ]lcavol X.scaled[train, ]lweight     X.scaled[train, ]age    X.scaled[train, ]lbph 
          2.46493292               0.67952814               0.26305307              -0.14146483               0.21014656 
X.scaled[train, ]svi     X.scaled[train, ]lcp X.scaled[train, ]gleason   X.scaled[train, ]pgg45 
          0.30520060              -0.28849277              -0.02130504               0.26695576 
| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ +1, my mistake was scaling y as well. $\endgroup$ – Gordon Gustafson Jan 12 at 22:16

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.