LOG for Stat 705 Course of Stepwise Model Fitting 4/7/08 ================================================= Start in R with "attitude" data, 30 x 7 matrix with columns > names(attitude) [1] "rating" "complaints" "privileges" "learning" "raises" [6] "critical" "advance" ## Use "rating" as response variable > tmpfit = lm( rating ~ . , data=attitude) > summary(tmpfit)$coef Estimate Std. Error t value Pr(>|t|) (Intercept) 10.78707639 11.5892572 0.9307824 0.3616337210 complaints 0.61318761 0.1609831 3.8090182 0.0009028679 privileges -0.07305014 0.1357247 -0.5382229 0.5955939205 learning 0.32033212 0.1685203 1.9008516 0.0699253459 raises 0.08173213 0.2214777 0.3690310 0.7154800884 critical 0.03838145 0.1469954 0.2611064 0.7963342642 advance -0.21705668 0.1782095 -1.2179862 0.2355770486 > extractAIC(tmpfit) [1] 7.0000 123.3635 ### first argument is degrees of freedom used > extractAIC(update(tmpfit, .~ complaints)) [1] 2.0000 118.6275 > summary(update(tmpfit, .~ complaints))$coef Estimate Std. Error t value Pr(>|t|) (Intercept) 14.3763194 6.6199860 2.171654 3.850833e-02 complaints 0.7546098 0.0975329 7.736978 1.987682e-08 extractAIC(update(tmpfit, .~ complaints)) ### Let's see what we would get "stepwise" using all initial ## variables plus pairwise interactions. > stpfit1 = step(tmpfit, ratings ~ (complaints + privileges + learning + raises + critical + advance)^2, direction = "forw") Start: AIC= 123.36 rating ~ complaints + privileges + learning + raises + critical + advance Df Sum of Sq RSS AIC + learning:advance 1 183.38 965.62 120.15 + complaints:advance 1 131.11 1017.89 121.73 1149.00 123.36 + learning:critical 1 72.80 1076.20 123.40 + privileges:advance 1 69.20 1079.80 123.50 + learning:raises 1 39.67 1109.34 124.31 + privileges:learning 1 17.63 1131.37 124.90 + complaints:learning 1 15.44 1133.56 124.96 + complaints:raises 1 14.26 1134.74 124.99 + critical:advance 1 13.50 1135.50 125.01 + raises:advance 1 13.03 1135.97 125.02 + complaints:critical 1 4.46 1144.54 125.25 + complaints:privileges 1 4.04 1144.96 125.26 + privileges:raises 1 3.56 1145.44 125.27 + privileges:critical 1 2.14 1146.86 125.31 + raises:critical 1 1.30 1147.70 125.33 Step: AIC= 120.15 rating ~ complaints + privileges + learning + raises + critical + advance + learning:advance Df Sum of Sq RSS AIC 965.62 120.15 + critical:advance 1 11.42 954.20 121.79 + complaints:privileges 1 9.98 955.64 121.84 + complaints:advance 1 7.94 957.69 121.90 + raises:advance 1 6.49 959.13 121.94 + raises:critical 1 5.35 960.27 121.98 + privileges:advance 1 1.40 964.22 122.10 + complaints:raises 1 1.06 964.56 122.11 + privileges:critical 1 1.06 964.56 122.11 + learning:critical 1 1.03 964.59 122.12 + complaints:learning 1 0.90 964.72 122.12 + complaints:critical 1 0.12 965.51 122.14 + privileges:raises 1 0.11 965.51 122.14 + privileges:learning 1 0.09 965.53 122.14 + learning:raises 1 0.04 965.59 122.15 > stpfit2 = step(lm(rating ~ 1, data=attitude), scope= ratings ~ (complaints + privileges + learning + raises + critical + advance)^2, direction = "forw") ... Step: AIC= 118 rating ~ complaints + learning Df Sum of Sq RSS AIC 1254.65 118.00 + advance 1 75.54 1179.11 118.14 + privileges 1 30.03 1224.62 119.28 + complaints:learning 1 22.55 1232.10 119.46 + raises 1 1.19 1253.46 119.97 + critical 1 0.002134 1254.65 120.00 > summary(lm(rating ~ complaints + learning, data=attitude))$coef Estimate Std. Error t value Pr(>|t|) (Intercept) 9.8708805 7.0612236 1.397899 1.735231e-01 complaints 0.6435176 0.1184774 5.431563 9.573795e-06 learning 0.2111918 0.1344037 1.571324 1.277539e-01 > stpfit3 = step(lm(rating ~ 1, data=attitude), scope= ratings ~ (complaints + privileges + learning + raises + critical + advance)^2, direction = "both", k=3.84) Step: AIC= 122.31 #### NOTE different criterion because different k rating ~ complaints Df Sum of Sq RSS AIC 1369.4 122.3 + learning 1 114.7 1254.6 123.5 + raises 1 11.1 1358.3 125.9 + privileges 1 7.5 1361.9 126.0 + advance 1 4.2 1365.2 126.1 + critical 1 0.009975 1369.4 126.1 - complaints 1 2927.6 4297.0 152.8 > write.table(attitude, "attitude.txt", row.names=F) ### delete "" using wordprocessor. ### Now try same analyses in SAS libname home "."; options nocenter linesize=80; data attitude; infile "attitude.txt"; input rating complaints privileges learning raises critical advance; if rating ne .; run; proc reg data = attitude; model rating = complaints privileges learning raises critical advance / selection = forward sle = .05; run; The REG Procedure Model: MODEL1 Dependent Variable: rating Forward Selection: Step 1 Parameter Standard Variable Estimate Error Type II SS F Value Pr > F Intercept 14.37632 6.61999 230.64710 4.72 0.0385 complaints 0.75461 0.09753 2927.58425 59.86 <.0001 Bounds on condition number: 1, 1 ----------------------------------------------------------------------- No other variable met the 0.0500 significance level for entry into the model. proc corr data=attitude; var rating; with privileges learning raises critical advance ; partial complaints; run; The CORR Procedure Pearson Partial Correlation Coefficients, N = 30 Prob > |r| under H0: Partial Rho=0 rating privileges -0.07410 0.7025 learning 0.28946 0.1278 raises 0.09004 0.6423 critical 0.00270 0.9889 advance -0.05505 0.7767