Module 11 assignment

 When I ran the different model formulas, I noticed that they all use the same factors, but depending on how I write the model, R builds the design matrix differently. So a*b gives main effects and the interaction. a:b only gives the interaction cells. a + b only gives main effects. Same variables, but different structure.

The main thing I saw is that some of these versions actually run fine, and some end up with singularities. For example, lm(z ~ a*b) worked normally, but lm(z ~ a:b) gave me an NA because the interaction only model wasn’t full rank with the intercept. So it isn’t the data that causes NAs, it’s literally how the model is written.

The implication is just that the way you write the formula matters. It controls the math in the background, and that decides whether R can estimate all the coefficients or not.

> a <- gl(2,2,8)   # same thing as g1(2,2,8)
> b <- gl(2,4,8)   # same thing as g1(2,4,8)
> x <- 1:8
> y <- c(1:4,8:5)
> set.seed(1)
> z <- rnorm(8)
> > model.matrix(~ a*b)
  (Intercept) a2 b2 a2:b2
1           1  0  0     0
2           1  0  0     0
3           1  1  0     0
4           1  1  0     0
5           1  0  1     0
6           1  0  1     0
7           1  1  1     1
8           1  1  1     1
attr(,"assign")
[1] 0 1 2 3
attr(,"contrasts")
attr(,"contrasts")$a
[1] "contr.treatment"

attr(,"contrasts")$b
[1] "contr.treatment"

> model.matrix(~ a:b)
  (Intercept) a1:b1 a2:b1 a1:b2 a2:b2
1           1     1     0     0     0
2           1     1     0     0     0
3           1     0     1     0     0
4           1     0     1     0     0
5           1     0     0     1     0
6           1     0     0     1     0
7           1     0     0     0     1
8           1     0     0     0     1
attr(,"assign")
[1] 0 1 1 1 1
attr(,"contrasts")
attr(,"contrasts")$a
[1] "contr.treatment"

attr(,"contrasts")$b
[1] "contr.treatment"

> model.matrix(~ a + b)
  (Intercept) a2 b2
1           1  0  0
2           1  0  0
3           1  1  0
4           1  1  0
5           1  0  1
6           1  0  1
7           1  1  1
8           1  1  1
attr(,"assign")
[1] 0 1 2
attr(,"contrasts")
attr(,"contrasts")$a
[1] "contr.treatment"

attr(,"contrasts")$b
[1] "contr.treatment"

> model.matrix(~ a*b - 1)  # no intercept
  a1 a2 b2 a2:b2
1  1  0  0     0
2  1  0  0     0
3  0  1  0     0
4  0  1  0     0
5  1  0  1     0
6  1  0  1     0
7  0  1  1     1
8  0  1  1     1
attr(,"assign")
[1] 1 1 2 3
attr(,"contrasts")
attr(,"contrasts")$a
[1] "contr.treatment"

attr(,"contrasts")$b
[1] "contr.treatment"

> > lm(z ~ a*b)

Call:
lm(formula = z ~ a * b)

Coefficients:
(Intercept)           a2           b2  
   -0.22141      0.60123     -0.02408  
      a2:b2  
    0.25713  

> lm(z ~ a:b)

Call:
lm(formula = z ~ a:b)

Coefficients:
(Intercept)        a1:b1        a2:b1  
     0.6129      -0.8343      -0.2331  
      a1:b2        a2:b2  
    -0.8584           NA  
Call:
lm(formula = vas ~ subject + period + treat, data = ash.long)

Residuals:
   Min     1Q Median     3Q    Max 
-48.94 -18.44   0.00  18.44  48.94 

Coefficients: (1 not defined because of singularities)
            Estimate Std. Error t value
(Intercept)  -113.06      27.39  -4.128
subject2       51.50      37.58   1.370
subject3      121.50      37.58   3.233
subject4       97.00      37.58   2.581
subject5      125.00      37.58   3.326
subject6       31.50      37.58   0.838
subject7      119.50      37.58   3.180
subject8      132.00      37.58   3.513
subject9       80.50      37.58   2.142
subject10     116.00      37.58   3.087
subject11     121.50      37.58   3.233
subject12     154.50      37.58   4.111
subject13     131.00      37.58   3.486
subject14     125.00      37.58   3.326
subject15      99.00      37.58   2.634
subject16      80.50      37.58   2.142
period2           NA         NA      NA
treat         -42.87      13.29  -3.227
            Pr(>|t|)    
(Intercept) 0.000895 ***
subject2    0.190721    
subject3    0.005573 ** 
subject4    0.020867 *  
subject5    0.004604 ** 
subject6    0.415070    
subject7    0.006215 ** 
subject8    0.003142 ** 
subject9    0.049003 *  
subject10   0.007518 ** 
subject11   0.005573 ** 
subject12   0.000925 ***
subject13   0.003318 ** 
subject14   0.004604 ** 
subject15   0.018768 *  
subject16   0.049003 *  
period2           NA    
treat       0.005644 ** 
---
Signif. codes:  
  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1
  ‘ ’ 1

Residual standard error: 37.58 on 15 degrees of freedom
Multiple R-squared:  0.7566,	Adjusted R-squared:  0.4969 
F-statistic: 2.914 on 16 and 15 DF,  p-value: 0.02229

>

Comments

Post a Comment

Popular posts from this blog

Module 10

 lm(cystfibr$spemax ~ age + weight + bmp + fev1, data=cystfibr) anova(lm(cystfibr$spemax ~ age + weight + bmp + fev1, data=cystfibr))   I used the cystfibr dataset from the ISwR package In completing this assignment. I ran a regression on spemax Based on the predictor variables age, weight, BMI an FEV 1 with spemex is being predicted. My results indicated that my overall model model was statistically significant. Therefore, at least one of the predictor variables is significantly associated with spemax .   Based upon the coefficients produced in this analysis, I examined age. Since the coefficient for age was negative, it can be concluded that the age increases spemex We'll decrease while controlling for the effects of the other predictor variables. This indicates that older patients generally exhibit lower maximum specific power.  
Read more

basic graph

Image
 For this assignment I created a histogram in R using the built in mtcars dataset to show the distribution of miles per gallon (mpg). The visualization makes it easy to see that most cars fall between about 15 and 25 mpg, with only a few vehicles above 30 mpg. This fits with Few’s and Yau’s discussion on basic descriptive visualization because it focuses on clarity and simplicity rather than decoration. The chart communicates one clear message about the spread of the data and avoids unnecessary visual elements. Overall, it works well as a basic example of descriptive analysis in R.
Read more

Module 3 assignemnt

  1. Central Tendency Set #1: 10, 2, 3, 2, 4, 2, 5 Mean = 4.0 Median = 3.0 Mode = 2 Set #2: 20, 12, 13, 12, 14, 12, 15 Mean = 14.0 Median = 13.0 Mode = 12 2. Variation Set #1 Range = 8 IQR = 2.5 Variance = 8.33 Standard Deviation = 2.89 Set #2 Range = 8 IQR = 2.5 Variance = 8.33 Standard Deviation = 2.89 3. Comparison Set #1 stays pretty low overall with an average of 4, a median of 3, and a mode of 2. Set #2 is higher across the board with an average around 14, a median of 13, and a mode of 12. Basically, Set #2 looks like Set #1 just shifted up by about 10 points. The spread is the same in both sets though. They have the exact same range, IQR, variance, and standard deviation. That means the consistency of the data is identical, even though one set has smaller values and the other has larger ones. So overall, the difference is just where the numbers sit . Set #1 clusters around lower values, while Set #2 clusters higher...
Read more