###Hierarchical regression, Parkinson's disease example

#First build a data set with just the needed variables and omit missing cases

Wk11a_nomiss<-na.omit(Wk11a[,c(75,14,98,27,28)])

#create two linear models 
#one block model
oneBlockModel <- lm(CC_PDQ_mob ~ PD_hoehn + mike_age, Wk11a_nomiss)
#two block model
twoBlockModel <- lm(CC_PDQ_mob ~ PD_hoehn + mike_age + C_falls + FinalTUG, Wk11a_nomiss)

#get summary data for each model using summary(OBJECT)
summary(oneBlockModel)

## 
## Call:
## lm(formula = CC_PDQ_mob ~ PD_hoehn + mike_age, data = Wk11a_nomiss)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -41.456  -7.045  -2.066   6.504  31.814 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -6.19313    1.44375  -4.290 1.88e-05 ***
## PD_hoehn     4.55000    0.21602  21.063  < 2e-16 ***
## mike_age     0.10998    0.02167   5.074 4.25e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 9.511 on 1997 degrees of freedom
## Multiple R-squared:  0.2129, Adjusted R-squared:  0.2121 
## F-statistic: 270.1 on 2 and 1997 DF,  p-value: < 2.2e-16

summary(twoBlockModel)

## 
## Call:
## lm(formula = CC_PDQ_mob ~ PD_hoehn + mike_age + C_falls + FinalTUG, 
##     data = Wk11a_nomiss)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -23.344  -5.541  -1.640   5.035  35.462 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -2.15171    1.21463  -1.771   0.0766 .  
## PD_hoehn     1.82200    0.20279   8.985   <2e-16 ***
## mike_age    -0.02566    0.01906  -1.346   0.1783    
## C_falls      2.61387    0.19931  13.115   <2e-16 ***
## FinalTUG     0.71739    0.03168  22.644   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 7.939 on 1995 degrees of freedom
## Multiple R-squared:  0.4522, Adjusted R-squared:  0.4511 
## F-statistic: 411.6 on 4 and 1995 DF,  p-value: < 2.2e-16

varoneblock<-summary(oneBlockModel)$r.squared
vartwoblock<-summary(twoBlockModel)$r.squared

#compare r-squareds between models (compute r-squared change)

(rsq.change1<-vartwoblock - varoneblock )

## [1] 0.2392263

#compare successive models using anova(MODEL1, MODEL2)
anova(oneBlockModel, twoBlockModel)

## Analysis of Variance Table
## 
## Model 1: CC_PDQ_mob ~ PD_hoehn + mike_age
## Model 2: CC_PDQ_mob ~ PD_hoehn + mike_age + C_falls + FinalTUG
##   Res.Df    RSS Df Sum of Sq      F    Pr(>F)    
## 1   1997 180659                                  
## 2   1995 125748  2     54911 435.58 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1