0
$\begingroup$

I have got 2 models. The first model uses a factor as the independent variable. The 2nd model uses a matrix with 0/1 coefficients as the independent variable.

# fit a model where an iris' petal width depends on species
# first up with a categorical variable
summary(lm(Petal.Width ~ 0+Species, data=iris))

# 2nd with a matrix indicating the species
smat = model.matrix(~ 0+Species, data=iris)
summary(lm(Petal.Width ~ 0+smat, data=iris))

The outcome of the two models above is the same. This is what I would expect:

Call:
lm(formula = Petal.Width ~ 0 + smat, data = iris)

Residuals:
   Min     1Q Median     3Q    Max 
-0.626 -0.126 -0.026  0.154  0.474 

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
smatSpeciessetosa      0.24600    0.02894    8.50 1.96e-14 ***
smatSpeciesversicolor  1.32600    0.02894   45.82  < 2e-16 ***
smatSpeciesvirginica   2.02600    0.02894   70.00  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.2047 on 147 degrees of freedom
Multiple R-squared:  0.9796,    Adjusted R-squared:  0.9792 
F-statistic:  2357 on 3 and 147 DF,  p-value: < 2.2e-16

However, when the model gets more complicated the summary is not the same anymore:

x = c(rep("x1",19),rep("x2",23),rep("x3",8),rep("x4",101),rep("x5",66),rep("x6",26),rep("x7",88))
y = c(0.74, 0.37, -0.34, 0.18, 0.18, -1.2, 0.04, 0.61, 0.85, 0, -0.03,  -0.18, -0.12, 1.37, -0.34, -0.4, -1.71, -0.13, 3.63, 1.6, 0.47,  0.68, 0.14, 0.86, 0.68, 0.07, -0.33, 4.14, 5.3, 0.94, -0.54,  -2.14, 0.12, 0.81, 1.51, -1.53, 3.37, -1.83, -6.99, -4.69, 1.59,  -1.01, -0.09, 0.04, 0.01, 20.35, -0.28, -1.67, -1.08, -5.33,  1.26, 4.35, 0.07, 1.55, 0.65, 0.15, -0.08, -0.31, -1.61, 0.57,  -0.36, 0.47, 0.72, -0.66, 0.92, 0.77, -1.93, -1.06, -0.93, 0.33,  0, 1.23, -1.53, -0.15, -3.5, 1.01, 0.63, 0.98, 1.06, 0.77, 0.33,  0.72, 1.56, 1.81, 0.74, 1.06, 0.09, 0.15, 0.75, 0.17, 0.77, 0.63,  0.54, 0.22, 0.49, -0.81, -0.03, 2.19, 0.76, 0.11, 0.94, -1.7,  -2.11, 0.32, -4.9, 0.39, 1.59, -5.66, 2.08, -10.43, -0.08, -0.9,  -0.16, -0.86, 0.43, -2.37, -2.16, -0.5, -0.6, 0.27, -0.34, 1.41,  -0.42, -1.43, -2.25, -0.69, -0.71, -1.54, -3.12, -0.93, -2.01,  0.26, -0.53, -1.09, -1.45, -2.75, -1.63, 0.16, -0.13, 0.08, -2.93,  -1.46, 0.18, -0.45, -0.88, 0.05, 0.68, -0.62, -2.55, -3.33, -2.35,  -3.56, -0.07, -0.29, 1.81, 1.07, 0.01, 0.85, -0.21, -1.68, 1.13,  1.58, -0.31, 1.14, -0.73, 2.04, 0.82, 0.78, 0.69, 1.18, 3.78,  0.97, 0.7, 13.63, -2.35, 4.52, 7.46, 4.09, 1.36, 0.77, -0.51,  0.09, 3.12, -3.36, -1.04, 2.79, 0.1, 2.83, 1.96, 0.15, -0.09,  0.39, 0.65, -1.7, -0.84, 1.86, -0.6, 0.18, 0.49, -0.66, -0.16,  2.16, 0.46, 1.66, 0.13, 2.12, -0.03, -0.31, 0.13, -0.45, 0.18,  -0.17, 0.25, -0.16, -0.06, 0.01, 0.83, -2.04, 1.06, 0.35, 0.86,  0.8, 2.01, 2.55, 0.07, 0.55, 1.59, -2.79, -2.77, 1.15, 0.46,  0.45, -0.4, 1.54, 0.29, 0.17, -0.36, 0.99, -9.59, 2.6, 1.79,  -0.56, -0.49, 0.28, 1.46, -1.19, -0.48, -2.52, 1.08, -4.1, -2.91,  -4.89, -0.12, -10.36, -4.97, -7.07, -3.81, -4.34, -12.82, -2.51,  -8.47, -6.12, -1.56, -4.48, -2.88, -1.94, -6.01, -1.35, 1.23,  3.08, -2.3, -0.99, 0.63, 2.4, 2.27, -0.79, 1.76, 1.21, -4.43,  -1.34, -1.94, -0.57, -2.77, -3.07, -1.15, -2.28, 0.68, 1.5, -1.78,  0.46, 0.39, -2.04, -1.22, 0.66, 1.22, 2.61, 5.21, 1.17, 2.73,  2.48, -0.95, 1.59, 0.14, -2.62, -1.73, -0.55, -1.25, 0.15, -0.11,  -0.23, -3.47, 1.46, -1.32, 0.24, 1.48, 2.68, 0.07, 0.24, -0.76,  -1.18, -1.57, -4.08, -2.74, 0.03, 0.34, 0.34, 0.84, 1.26, -0.81,  -0.79, -1.59)

with the following output:

> summary(lm(y ~ 0+x))

Call:
lm(formula = y ~ 0 + x)

Residuals:
     Min       1Q   Median       3Q      Max 
-11.6615  -1.0202   0.0515   1.1498  18.8562 

Coefficients:
    Estimate Std. Error t value Pr(>|t|)    
xx1  0.18526    0.59168   0.313    0.754    
xx2  0.14000    0.53778   0.260    0.795    
xx3  1.49375    0.91184   1.638    0.102    
xx4 -0.41149    0.25663  -1.603    0.110    
xx5  0.81182    0.31746   2.557    0.011 *  
xx6  0.01077    0.50580   0.021    0.983    
xx7 -1.15852    0.27493  -4.214 3.26e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.579 on 324 degrees of freedom
Multiple R-squared:  0.08401,   Adjusted R-squared:  0.06422 
F-statistic: 4.245 on 7 and 324 DF,  p-value: 0.0001676

> summary(lm(y ~ 0+model.matrix(~x)))

Call:
lm(formula = y ~ 0 + model.matrix(~x))

Residuals:
     Min       1Q   Median       3Q      Max 
-11.6615  -1.0202   0.0515   1.1498  18.8562 

Coefficients:
                            Estimate Std. Error t value Pr(>|t|)  
model.matrix(~x)(Intercept)  0.18526    0.59168   0.313   0.7544  
model.matrix(~x)xx2         -0.04526    0.79956  -0.057   0.9549  
model.matrix(~x)xx3          1.30849    1.08699   1.204   0.2296  
model.matrix(~x)xx4         -0.59675    0.64494  -0.925   0.3555  
model.matrix(~x)xx5          0.62656    0.67147   0.933   0.3515  
model.matrix(~x)xx6         -0.17449    0.77841  -0.224   0.8228  
model.matrix(~x)xx7         -1.34379    0.65244  -2.060   0.0402 *
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.579 on 324 degrees of freedom
Multiple R-squared:  0.08401,   Adjusted R-squared:  0.06422 
F-statistic: 4.245 on 7 and 324 DF,  p-value: 0.0001676

Why are the p-values of the coefficients different here?

If it's the order in which the coefficients are fitted, is there a way to get the same behaviour?

Improve this question
$\endgroup$
4
  • 1
    $\begingroup$ You have, whether you knew it or not, used a different coding scheme for the predictors in the two models, so the coefficients have different interpretations. Also, please try to write questions that don't rely on everyone understanding R code. $\endgroup$
    – Scortchi ♦
    Commented May 13, 2015 at 12:02
  • $\begingroup$ I just realized that I made a mistake in the encoding of the coefficients in the last example - will see how this changes the outcome. $\endgroup$
    – Michael Schubert
    Commented May 13, 2015 at 12:10
  • 1
    $\begingroup$ Ok, so the answer is: I forgot that the command model.matrix(~x) will fit coefficients with intercept, which I didn't want. The right syntax is model.matrix(~0+x), and then the results are the same. (for R: I added 2 lines to explain the model fit) $\endgroup$
    – Michael Schubert
    Commented May 13, 2015 at 12:18
  • $\begingroup$ You have the ability to answer your own question. Doing so will show that this question has an answer, and improve the website for future visitors. $\endgroup$
    – Sycorax ♦
    Commented May 13, 2015 at 12:20

1 Answer 1

Reset to default
1
$\begingroup$

I forgot that the command

model.matrix(~x)

will fit coefficients with intercept, which I didn't want. The right syntax is 

model.matrix(~0+x)

and then the results are the same.

So, to sum up, the following 2 models behave the same way:

summary(lm(y ~ 0+x))
summary(lm(y ~ 0+model.matrix(~0+x))) # notice that '0' appears twice

To get the same behaviour with intercept, the following syntax is needed:

summary(lm(y ~ x))
summary(lm(y ~ model.matrix(~x))) # notice that no '0' is in model.matrix()
Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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