I came across a puzzling difference in standardized (beta) coefficients with linear regression model computed with R and SPSS using dummy coded variables. I have used the hsb2 data set and created a contrast (dummy coding), so that the third category is the reference. Here is the R code:
# Read the data
hsb2 <- read.table('https://stats.idre.ucla.edu/stat/data/hsb2.csv', header = TRUE, sep = ",")
# Create a factor variable with respondents' race
hsb2$race.f <- factor(hsb2$race, labels = c("Hispanic", "Asian", "African-Am", "Caucasian"))
# Add a contrast (dummy coding) to the new race variable, so that the third category is the reference.
contrasts(hsb2$race.f) <- contr.treatment(n = 4, base = 3)
# Scale the writing achievement score (mean of 0 and SD of 1), it will be the dependent variable
hsb2$write <- scale(hsb2$write)
# Fit the model and print the summary
summary(lm(write ~ race.f, hsb2))
The output I get:
Call:
lm(formula = write ~ race.f, data = hsb2)
Residuals:
Min 1Q Median 3Q Max
-2.43234300577889240 -0.57585945002954031 0.10259059641484436 0.73850677561040290 1.98341819735365221
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.48266692834536767 0.21290900103341129 -2.26700999999999997 0.0244812 *
race.f1 -0.18374751916973245 0.28828015018135283 -0.63739000000000001 0.5246133
race.f2 1.03390948585456388 0.35741973343705952 2.89270000000000005 0.0042513 **
race.f4 0.61772635713618673 0.22711822910747051 2.71984000000000004 0.0071181 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.050000000000000003 ‘.’ 0.10000000000000001 ‘ ’ 1
Residual standard error: 0.95215799866456285 on 196 degrees of freedom
Multiple R-squared: 0.1070625554447362515, Adjusted R-squared: 0.09339514557909434078
F-statistic: 7.833419535758452845 on 3 and 196 DF, p-value: 0.000057845156841983661
However, when I run the same analysis with SPSS I get quite different beta regression coefficients, here is the code:
* Create the dummy variables.
RECODE race (1 = 1) (ELSE = 0) INTO race.f1.
RECODE race (2 = 1) (ELSE = 0) INTO race.f2.
RECODE race (3 = 1) (ELSE = 0) INTO race.f3.
RECODE race (4 = 1) (ELSE = 0) INTO race.f4.
EXECUTE.
* Execute the analysis, so that the third category is the reference.
REGRESSION
/MISSING LISTWISE
/STATISTICS COEFF OUTS R ANOVA
/CRITERIA=PIN(.05) POUT(.10)
/NOORIGIN
/DEPENDENT write
/METHOD=ENTER race.f1 race.f2 race.f4.
Here is the SPSS output I get:
What really baffles me is that the everything else is the same (model statistics - R2, adjusted R2, degrees of freedom, F-statistic; and the t-values and p-values of the beta regression coefficients), but the standardized beta regression coefficients are not even close. If I run without standardization, the unstandardized regression coefficients and all other statistics match between R and SPSS.
Can anyone help with this? Am I missing something?
EDIT Following the source provided by aosmith (thanks once again), I did the dummy coding by hand, scaling the separate dummies:
hsb2 <- read.table('https://stats.idre.ucla.edu/stat/data/hsb2.csv', header = TRUE, sep = ",")
hsb2$write <- scale(hsb2$write)
hsb2$race.f1 <- scale(hsb2$race == 1)
hsb2$race.f2 <- scale(hsb2$race == 2)
hsb2$race.f3 <- scale(hsb2$race == 3)
hsb2$race.f4 <- scale(hsb2$race == 4)
summary(lm(write ~ race.f1 + race.f2 + race.f4, hsb2))
I got exactly the same results as in SPSS:
Call:
lm(formula = write ~ race.f1 + race.f2 + race.f4, data = hsb2)
Residuals:
Min 1Q Median 3Q Max
-2.4323430057788924 -0.5758594500295402 0.1025905964148444 0.7385067756104029 1.9834181973536520
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.000000000000000030665367318040625 0.067327737761672404315227424831392 0.00000000000000000 1.0000000
race.f1 -0.059860715422078700220787084163021 0.093915042280922900186368451613816 -0.63739000000000001 0.5246133
race.f2 0.236302452210854940783946176452446 0.081689123308428354675037041943142 2.89270000000000005 0.0042513 **
race.f4 0.276515793804944842726456499804044 0.101666015515960786452787090183847 2.71984000000000004 0.0071181 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.050000000000000003 ‘.’ 0.10000000000000001 ‘ ’ 1
Residual standard error: 0.95215799866456285 on 196 degrees of freedom
Multiple R-squared: 0.1070625554447362238, Adjusted R-squared: 0.09339514557909434078
F-statistic: 7.833419535758451957 on 3 and 196 DF, p-value: 0.000057845156841983668
However, using this approach in a custom function would not be quite handy. I wonder if there is a way to do it still using the contrasts function to assign the dummies.
contrastsinstead of making a dummies by hand which would not be very handy when constructing a function or using other contrast coding schemes... Hm... – panman