0
$\begingroup$

To my knowledge, mean-centering does affect the values of regression coefficient for variables involved with interaction terms. But it does reduce the standard errors of the coefficient estimates for variables involved with interaction terms.

But my simulation results shows weird shift in regression coefficients after fitting a linear model with mean-centered dependent variables with interaction terms. Here's the simulation step:

  1. Generate dependent variable X from a multivariate normal distribution

  2. Generate independent variable y described by a linear model of:

$$y = 15 +13x_1 + 0.5 x_2 + 5 x_3 + 42 x_4 -2 x^2 + 9x_2x_3 + \varepsilon$$

from X, with gaussian noise ($\varepsilon$),

  1. Mean-center X and y

  2. Add interaction terms to X

  3. Regress mean-centered X against mean-centered y

Here's the Python code snippet

import numpy as np
import statsmodels.api as sm
import matplotlib.pyplot as plt
from sklearn.preprocessing import StandardScaler
import pandas as pd

######################################## Collinear random data generation #########################################

def generate_collinear_data(cov, n=10, true_coefs=[0, 0], true_intercept=0, feature_means=[0, 0, 0], loc=0, scale=1, 
                            true_coefs_interaction=[0, 0]):

    # random generation of 3D gaussian collinear features.
    X = np.random.multivariate_normal(mean=feature_means, cov=cov, size=n)

    # generate gaussian white noise.
    gaussian_noise = np.random.normal(loc=loc, scale=scale, size=n)

    # make the outcome.
    y = true_intercept + gaussian_noise
    for i in range(len(true_coefs)):
        y += true_coefs[i] * X[:, i]

    y += true_coefs_interaction[0] * (X[:, 0] ** 2) + true_coefs_interaction[1] * X[:, 1] * X[:, 2]

    return X, y

# settings
m = 1000                                   # number of simulations
correlation = 0.01                          # degree of collinearity (HIGH)

kwargs = {
    'n': 10000,                             # sample size
    'true_coefs': [13, 0.5, 5, 42],         # linear regression coefficients, 2 features
    'true_intercept': 15,                   # y-intercept
    'feature_means': [-12, 14, 2, -20],     # means multivariate normal distribution. This is not important
    'loc': 0,                               # mean of gaussian noise
    'scale': 1,                             # standard deviation of gaussian noise
    'true_coefs_interaction': [-2, 9]       # linear regression coefficients of interaction terms
}

# high collinearity covariance matrix
cov = np.full((len(kwargs['true_coefs']), len(kwargs['true_coefs'])), correlation)
np.fill_diagonal(cov, 1)

print('True intercept + coefficients:  ', [kwargs['true_intercept']] + kwargs['true_coefs'] + kwargs['true_coefs_interaction'])
print()

############################################# Original #############################################

X_org, y_org = generate_collinear_data(cov, **kwargs)

interaction_org_1 = X_org[:, 0].reshape(-1, 1) ** 2                            # x1^2
interaction_org_2 = X_org[:, 1].reshape(-1, 1) * X_org[:, 2].reshape(-1, 1)    # x2 * x3
X_org = np.concatenate([X_org, interaction_org_1, interaction_org_2], axis=1)

X_st_org = sm.add_constant(X_org)
model_org = sm.OLS(y_org, X_st_org).fit()

print('Original: ', [ '%.2f' % elem for elem in model_org.params ])

########################################## Mean-centered ###########################################

X_std, y_std = generate_collinear_data(cov, **kwargs)

# centering
scaler = StandardScaler(with_mean=True, with_std=False)
X_std, y_std = scaler.fit_transform(X_std), scaler.fit_transform(y_std.reshape(-1, 1)).flatten()

# add interaction terms
interaction_std_1 = X_std[:, 0].reshape(-1, 1) ** 2                            # x1^2
interaction_std_2 = X_std[:, 1].reshape(-1, 1) * X_std[:, 2].reshape(-1, 1)    # x2 * x3
X_std = np.concatenate([X_std, interaction_std_1, interaction_std_2], axis=1)

# fit
X_std = sm.add_constant(X_std)
model_std = sm.OLS(y_std, X_std).fit()

print('Centering:  ', [ '%.2f' % elem for elem in model_std.params ])

I ran two regressions: with centering, and without centering. Here's the output:

>>> True intercept + coefficients:   [15, 13, 0.5, 5, 42, -2, 9]

>>> Original:  ['14.87', '13.07', '0.53', '5.11', '41.99', '-2.00', '8.99']
>>> Centering:   ['1.87', '61.02', '18.47', '130.99', '42.02', '-2.00', '8.98']

I see that the regression coefficients for $x_1$, $x_2$, and $x_3$ with standardized dependent variables are very different from the true regression coefficients. This goes against the theory that "mean-centering does not affect the true estimation of the regression coefficient".

What's more interesting is, that my simulation results correctly guesses value of coefficients for interaction terms, but fails to guess the values of coefficients for the individual independent variables involved with the interaction terms.

Why is this happening?

$\endgroup$
  • 2
    $\begingroup$ I think this answer basically covers your issues. $\endgroup$ – Huy Pham Jan 13 at 8:37