> 
> n <- length(year)
> 
> casino.mod1 <- lm(profit ~ slot + table)
> summary(casino.mod1)

Call:
lm(formula = profit ~ slot + table)

Residuals:
    Min      1Q  Median      3Q     Max 
-3.6215 -0.1834  0.4233  0.8339  1.5642 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  0.19320    0.81761   0.236  0.81471    
slot         0.15902    0.02282   6.968  6.8e-08 ***
table        0.55179    0.16534   3.337  0.00216 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 1.307 on 32 degrees of freedom
Multiple R-squared: 0.8621,     Adjusted R-squared: 0.8534 
F-statistic:   100 on 2 and 32 DF,  p-value: 1.716e-14 

> anova(casino.mod1)
Analysis of Variance Table

Response: profit
          Df Sum Sq Mean Sq F value    Pr(>F)    
slot       1 322.61  322.61 188.863 5.679e-15 ***
table      1  19.02   19.02  11.137  0.002155 ** 
Residuals 32  54.66    1.71                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
> SSE.mod1 <- deviance(casino.mod1)
> e.mod1 <- residuals(casino.mod1)
> 
> plot(e.mod1,type="o",xlab="Time",ylab="Residual")
> 
> #### Compute the Durbin-Watson Stat and Cochrane-Orcutt Autocorrelation estimate
> 
> DW1.mod1 <- 0
> 
> CO1.mod1 <- 0
> CO2.mod1 <- 0
> 
> for (t in 2:n) {
+ DW1.mod1 <- DW1.mod1 + (e.mod1[t] - e.mod1[t-1])^2
+ CO1.mod1 <- CO1.mod1 + e.mod1[t-1]*e.mod1[t]
+ CO2.mod1 <- CO2.mod1 + (e.mod1[t-1])^2
+ }
> 
> (DW.mod1 <- DW1.mod1/SSE.mod1)
        2 
0.3442027 
> (COr.mod1 <- CO1.mod1/CO2.mod1)
       1 
0.929251 
> 
> ################### Cochrane-Orcutt Method
> 
> ###### Transformed profit, slot, table
> 
> Yt <- profit[2:n] - (COr.mod1 * profit[1:(n-1)])
> X1t <- slot[2:n] - (COr.mod1 * slot[1:(n-1)])
> X2t <- table[2:n] - (COr.mod1 * table[1:(n-1)])
> 
> 
> 
> CO.mod1t <- lm(Yt ~ X1t + X2t)
> (CO.mod1t.sum <- summary(CO.mod1t))

Call:
lm(formula = Yt ~ X1t + X2t)

Residuals:
    Min      1Q  Median      3Q     Max 
-2.3558 -0.3326  0.1087  0.4173  1.0425 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.53291    0.20857  -2.555   0.0157 *  
X1t          0.34148    0.05618   6.078 9.82e-07 ***
X2t          0.55161    0.24912   2.214   0.0343 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.6551 on 31 degrees of freedom
Multiple R-squared: 0.6309,     Adjusted R-squared: 0.6071 
F-statistic:  26.5 on 2 and 31 DF,  p-value: 1.950e-07 

> anova(CO.mod1t)
Analysis of Variance Table

Response: Yt
          Df  Sum Sq Mean Sq F value    Pr(>F)    
X1t        1 20.6371 20.6371  48.093 8.864e-08 ***
X2t        1  2.1039  2.1039   4.903   0.03430 *  
Residuals 31 13.3024  0.4291                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
> e.CO1t <- residuals(CO.mod1t)
> SSE.CO1t <- deviance(CO.mod1t)
> 
> ###### Compute the D-W Stat for transformed model
> 
> DW1.CO1t <- 0
> 
> for (t in 2:(n-1)) DW1.CO1t <- DW1.CO1t + (e.CO1t[t] - e.CO1t[t-1])^2
> (DW.CO1t <- DW1.CO1t/SSE.CO1t)
       2 
2.350875 
> 
> ###### Back-Transform the Coefficients
> 
> (b0.CO1o <- coef(CO.mod1t.sum)[1,1]/(1-COr.mod1))
       1 
-7.53244 
> (b1.CO1o <- coef(CO.mod1t.sum)[2,1])
[1] 0.3414778
> (b2.CO1o <- coef(CO.mod1t.sum)[3,1])
[1] 0.5516122
> 
> (s.b0.CO1o <- coef(CO.mod1t.sum)[1,2]/(1-COr.mod1))
       1 
2.948052 
> (s.b1.CO1o <- coef(CO.mod1t.sum)[2,2])
[1] 0.05617966
> (s.b2.CO1o <- coef(CO.mod1t.sum)[3,2])
[1] 0.2491170
> 
> 
> ##################### Hildreth-Lu Method #########
> 
> SSE.r <- matrix(rep(0,2*length(seq(0.01,0.99,0.01))),ncol=2)
> r.count <- 0
> HL.r <- 0
> min.SSE <- 10000000
> 
> for (r in seq(0.01,0.99,0.01))   {
+ r.count <- r.count + 1
+ 
+ Yt <- profit[2:n] - (r * profit[1:(n-1)])
+ X1t <- slot[2:n] - (r * slot[1:(n-1)])
+ X2t <- table[2:n] - (r * table[1:(n-1)])
+ 
+ SSE.r[r.count,1] <- r
+ SSE.r[r.count,2] <- deviance(lm(Yt ~ X1t + X2t))
+ if (deviance(lm(Yt ~ X1t + X2t)) < min.SSE) {
+ min.SSE <- deviance(lm(Yt ~ X1t + X2t))
+ HL.r <- r 
+ }
+ }
> 
> HL.r
[1] 0.97
> SSE.r
       [,1]     [,2]
  [1,] 0.01 53.61723
  [2,] 0.02 52.85343
  [3,] 0.03 52.09780
  [4,] 0.04 51.35033
  [5,] 0.05 50.61101
  [6,] 0.06 49.87985
  [7,] 0.07 49.15684
  [8,] 0.08 48.44198
  [9,] 0.09 47.73527
 [10,] 0.10 47.03670
 [11,] 0.11 46.34627
 [12,] 0.12 45.66399
 [13,] 0.13 44.98983
 [14,] 0.14 44.32381
 [15,] 0.15 43.66591
 [16,] 0.16 43.01614
 [17,] 0.17 42.37448
 [18,] 0.18 41.74093
 [19,] 0.19 41.11550
 [20,] 0.20 40.49816
 [21,] 0.21 39.88892
 [22,] 0.22 39.28777
 [23,] 0.23 38.69470
 [24,] 0.24 38.10971
 [25,] 0.25 37.53278
 [26,] 0.26 36.96391
 [27,] 0.27 36.40309
 [28,] 0.28 35.85030
 [29,] 0.29 35.30555
 [30,] 0.30 34.76881
 [31,] 0.31 34.24007
 [32,] 0.32 33.71932
 [33,] 0.33 33.20654
 [34,] 0.34 32.70172
 [35,] 0.35 32.20484
 [36,] 0.36 31.71588
 [37,] 0.37 31.23482
 [38,] 0.38 30.76164
 [39,] 0.39 30.29631
 [40,] 0.40 29.83880
 [41,] 0.41 29.38910
 [42,] 0.42 28.94716
 [43,] 0.43 28.51295
 [44,] 0.44 28.08643
 [45,] 0.45 27.66757
 [46,] 0.46 27.25631
 [47,] 0.47 26.85261
 [48,] 0.48 26.45641
 [49,] 0.49 26.06765
 [50,] 0.50 25.68627
 [51,] 0.51 25.31220
 [52,] 0.52 24.94535
 [53,] 0.53 24.58564
 [54,] 0.54 24.23297
 [55,] 0.55 23.88723
 [56,] 0.56 23.54830
 [57,] 0.57 23.21606
 [58,] 0.58 22.89036
 [59,] 0.59 22.57103
 [60,] 0.60 22.25790
 [61,] 0.61 21.95078
 [62,] 0.62 21.64943
 [63,] 0.63 21.35362
 [64,] 0.64 21.06307
 [65,] 0.65 20.77750
 [66,] 0.66 20.49655
 [67,] 0.67 20.21987
 [68,] 0.68 19.94704
 [69,] 0.69 19.67763
 [70,] 0.70 19.41112
 [71,] 0.71 19.14700
 [72,] 0.72 18.88466
 [73,] 0.73 18.62348
 [74,] 0.74 18.36278
 [75,] 0.75 18.10185
 [76,] 0.76 17.83995
 [77,] 0.77 17.57632
 [78,] 0.78 17.31023
 [79,] 0.79 17.04100
 [80,] 0.80 16.76803
 [81,] 0.81 16.49089
 [82,] 0.82 16.20938
 [83,] 0.83 15.92363
 [84,] 0.84 15.63420
 [85,] 0.85 15.34221
 [86,] 0.86 15.04944
 [87,] 0.87 14.75845
 [88,] 0.88 14.47263
 [89,] 0.89 14.19620
 [90,] 0.90 13.93410
 [91,] 0.91 13.69186
 [92,] 0.92 13.47519
 [93,] 0.93 13.28969
 [94,] 0.94 13.14033
 [95,] 0.95 13.03103
 [96,] 0.96 12.96434
 [97,] 0.97 12.94124
 [98,] 0.98 12.96109
 [99,] 0.99 13.02182
> 
> Yt <- profit[2:n] - (HL.r * profit[1:(n-1)])
> X1t <- slot[2:n] - (HL.r * slot[1:(n-1)])
> X2t <- table[2:n] - (HL.r * table[1:(n-1)])
> 
> 
> HL.mod1t <- lm(Yt ~ X1t + X2t)
> (HL.mod1t.sum <- summary(HL.mod1t))

Call:
lm(formula = Yt ~ X1t + X2t)

Residuals:
    Min      1Q  Median      3Q     Max 
-2.3402 -0.3436  0.0385  0.4463  0.9578 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.34724    0.14435  -2.406   0.0223 *  
X1t          0.35681    0.05691   6.270 5.71e-07 ***
X2t          0.42980    0.23145   1.857   0.0728 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.6461 on 31 degrees of freedom
Multiple R-squared: 0.6529,     Adjusted R-squared: 0.6306 
F-statistic: 29.16 on 2 and 31 DF,  p-value: 7.518e-08 

> anova(HL.mod1t)
Analysis of Variance Table

Response: Yt
          Df  Sum Sq Mean Sq F value    Pr(>F)    
X1t        1 22.9085 22.9085 54.8759 2.421e-08 ***
X2t        1  1.4395  1.4395  3.4483   0.07284 .  
Residuals 31 12.9412  0.4175                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
> e.HL1t <- residuals(HL.mod1t)
> SSE.HL1t <- deviance(HL.mod1t)
> 
> ##### Compute D-W Stat
> 
> DW1.HL1t <- 0
> 
> for (t in 2:(n-1)) DW1.HL1t <- DW1.HL1t + (e.HL1t[t] - e.HL1t[t-1])^2
> (DW.HL1t <- DW1.HL1t/SSE.HL1t)
       2 
2.451622 
> 
> ##### Back-transform the coefficients
> 
> (b0.HL1o <- coef(HL.mod1t.sum)[1,1]/(1-HL.r))
[1] -11.5748
> (b1.HL1o <- coef(HL.mod1t.sum)[2,1])
[1] 0.3568061
> (b2.HL1o <- coef(HL.mod1t.sum)[3,1])
[1] 0.4297994
> 
> (s.b0.HL1o <- coef(HL.mod1t.sum)[1,2]/(1-HL.r))
[1] 4.811782
> (s.b1.HL1o <- coef(HL.mod1t.sum)[2,2])
[1] 0.05690984
> (s.b2.HL1o <- coef(HL.mod1t.sum)[3,2])
[1] 0.2314517
> 
> 
> ############## First Differences Method
> 
> Yt <- profit[2:n]-profit[1:(n-1)]
> X1t <- slot[2:n]-slot[1:(n-1)]
> X2t <- table[2:n]-table[1:(n-1)]
> 
> FD.mod1t <- lm(Yt ~ -1 + X1t + X2t)
> (FD.mod1t.sum <- summary(FD.mod1t))

Call:
lm(formula = Yt ~ -1 + X1t + X2t)

Residuals:
    Min      1Q  Median      3Q     Max 
-2.4904 -0.5461 -0.1671  0.2699  0.7592 

Coefficients:
    Estimate Std. Error t value Pr(>|t|)    
X1t  0.33063    0.05714   5.786 2.02e-06 ***
X2t  0.28906    0.23006   1.256    0.218    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.6765 on 32 degrees of freedom
Multiple R-squared: 0.6267,     Adjusted R-squared: 0.6033 
F-statistic: 26.86 on 2 and 32 DF,  p-value: 1.425e-07 

> anova(FD.mod1t)
Analysis of Variance Table

Response: Yt
          Df  Sum Sq Mean Sq F value    Pr(>F)    
X1t        1 23.8598 23.8598 52.1339 3.351e-08 ***
X2t        1  0.7225  0.7225  1.5788    0.2180    
Residuals 32 14.6452  0.4577                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
> 
> (b0.FD1o <- mean(profit) - (coef(FD.mod1t.sum)[1,1]*mean(slot)) 
+    - (coef(FD.mod1t.sum)[2,1]*mean(table)))
[1] -2.151775
> 
> (b1.FD1o <- coef(FD.mod1t.sum)[1,1])
[1] 0.3306297
> (b2.FD1o <- coef(FD.mod1t.sum)[2,1])
[1] 0.2890628
> 
> (s.b1.FD1o <- coef(FD.mod1t.sum)[1,2])
[1] 0.05714351
> (s.b2.FD1o <- coef(FD.mod1t.sum)[2,2])
[1] 0.2300562
> 
> 
> 
>