Recovering original regression coefficients from standardized

Question

Suppose I use Least Squares to estimate coefficients in the standard linear model with design matrix $X$'s columns standardized, so the model is $$ E[y] = X^*\beta^* $$ where $X^*$ is $X$ with columns centered and scaled. Assume $X$ has full column rank. To recover the regression coefficients $\beta$ in the model with unstandardized X $$ E[y] = X\beta $$ I should be able to use the equation $$ \qquad \quad X\beta = X^*\beta^* \\ \implies X^T X\beta = X^TX^*\beta^* \\ \qquad \quad \; \; \; \; \implies \qquad \beta = (X^T X)^{-1}X^TX^*\beta^* $$ To test this I ran some simple R code:

set.seed(100)
# set number of samples
N <- 10
# set number of regressors
p <- 3
y <- runif(N)
X <- matrix(runif(N*p), N)
# Least squares coefficients with unstandardized X
(b1 <- solve(crossprod(X), crossprod(X,y)))

# now standardize X and get new coefficients
Xs <- scale(X)
(b2 <- solve(crossprod(Xs), crossprod(Xs,y)))

# b.orig should be exactly the same as b1, but its not!
(b.orig <- solve(crossprod(X), crossprod(X, Xs %*% b2)))

Here is the output of b1:

            [,1]
[1,]  0.17109189
[2,]  0.52204169
[3,] -0.02115178

and the output of b.orig, which should equal b1, but does not:

            [,1]
[1,] -0.15935376
[2,]  0.16915433
[3,] -0.04696604

What is going wrong here? I have already looked at:http://bit.ly/1HxHS6U, who seems to use a similar idea (equating expected or fitted values of y in both equations).

Answer 1

The models are not the same. Therefore the coefficients should differ. When you recenter, you are taking linear combinations of the columns of $X$ with the vector $\mathbf{1}=(1,1,\ldots, 1)^\prime$. This is fine, provided that $\mathbf{1}$ lies in the column space. In your example it does not.

What is worse, when you do include $\mathbf{1}$ as a column, the entire calculation falls apart due to singularities.

In detail, the original model (including a constant term) should be

$$\mathbb{E}(Y) = \beta_0 + \beta_1 X_1 + \cdots + \beta_p X_p = \beta_0 + \sum_{i=1}^p \beta_i X_i.$$

Standardizing to $X_i = \sigma_i Z_i + \mu_i$ yields

$$\mathbb{E}(Y) = \beta_0 + \sum_{i=1}^p \beta_i (\sigma_i Z_i + \mu_i) = \left(\beta_0 + \sum_{i=1}^p \beta_i \mu_i\right) + \sum_{i=1}^p (\beta_i \sigma_i) Z_i = \beta_0^{*} + \sum_{i=1}^p \beta_i^{*} Z_i,$$

with $\beta_i^{*} = \sigma_i \beta_i$ giving the correct relationships between the "standardized" and unstandardized coefficients. (The usual definition of standardized coefficient also standardizes the response variable, so these $\beta_i^{*}$ might better be called "semi-standardized.")