I have data that require scaling fixed effects to fit linear mixed models with the lmer function from the lme4 package that converge. I want to unscale the estimated coefficients to present results in useful units. Previous posts address this issue, but I am unable to successfully unscale coefficients from models including categorical predictors and interaction terms. I used lmer() with my data set, but I have created an example using lm and the diamonds data for simplicity.
related threads: https://stackoverflow.com/questions/23642111/how-to-unscale-the-coefficients-from-an-lmer-model-fitted-with-a-scaled-respon
and https://stackoverflow.com/questions/24268031/unscale-and-uncenter-glmer-parameters
Below is some code that runs without error, but demonstrates that I can't get unscaling to work correctly:
library(ggplot2) library(dplyr) # function to rescale variables around 0 zscore <- function(x) (x - mean(x)) / sd(x) # creates new dataframe with scaled values included diamonds_s <- diamonds %>% select(carat, cut, price, depth) %>% mutate_each(funs(s = zscore(.)), -cut) Create unscaled and scaled models
# Linear model with unscaled interaction specified m1 <- lm(price ~ carat * depth * cut, data= diamonds_s) b1<- coef(m1) # summary(m1) # Linear model with scaled interaction specified m1s <- lm(price_s ~ carat_s * depth_s * cut, data= diamonds_s) b1s <- coef(m1s) # summary(m1s) Function to rescale coefficients from answer by Ben Bolker
# https://stackoverflow.com/questions/24268031/unscale-and-uncenter-glmer-parameters rescale.coefs <- function(beta,mu,sigma) { beta2 <- beta ## inherit names etc. beta2[-1] <- sigma[1]*beta[-1]/sigma[-1] beta2[1] <- sigma[1]*beta[1]+mu[1]-sum(beta2[-1]*mu[-1]) beta2 } Use model.matrix to get predictor variables
X <- model.matrix(price ~ carat * depth * cut, data= diamonds_s) # matrix of response (m[1]) and predictor (m[-1]) means m <- c(mean(diamonds_s$price), colMeans(X)[-1] ) s <- c(sd(diamonds_s$price), apply(X, 2, sd)[-1] ) # or with unscaled response (not run because I scaled the response, price_s): # m <- c(0, colMeans(X)[-1] ) # same for standard deviations # s <- c(1, apply(X, 2, sd)[-1] ) Now attempt to unscale
b1r <- rescale.coefs(b1s, m, s) The values obtained are wrong, however.
all.equal(b1, b1r, check.names=FALSE, tolerance=0.001) I appreciate any input about this issue.
