-------------
      name:  <unnamed>
       log:  /Users/eric/Work/MB/ATCR/200/Labs/lab 9.smcl
  log type:  smcl
 opened on:  15 Oct 2018, 08:36:25

.
. /* Suppose you run a moving company and want to predict the number of hours
> of labor for a prospective job.  You have collected some data on labor hours,
> number of rooms in the pickup location, and the volume of the items to be
> moved.  We will use these data to build a regression model for labor hours.
>
> We will read in the data two ways, first using the input command: */
.
. clear

. input laborhrs rooms volume

      laborhrs      rooms     volume
  1. 24.00           3.5     545
  2. 13.50           2.0     400
  3. 26.25           2.5     562
  4. 25.00           3.0     540
  5. 9.00            1.0     220
  6. 20.00           3.0     344
  7. 22.00           3.5     569
  8. 11.25           2.0     340
  9. 50.00           5.0     900
 10. 12.00           1.5     285
 11. 38.75           5.0     865
 12. 40.00           4.5     831
 13. 19.50           3.0     344
 14. 18.00           2.5     360
 15. 28.00           4.0     750
 16. 27.00           3.5     650
 17. 21.00           3.0     415
 18. 15.00           2.5     275
 19. 25.00           3.0     557
 20. 45.00           5.5     1028
 21. 29.00           4.5     793
 22. 21.00           3.0     523
 23. 22.00           3.5     564
 24. 16.50           2.5     312
 25. 37.00           4.0     757
 26. 32.00           3.5     600
 27. 34.00           4.0     796
 28. 25.00           3.5     577
 29. 31.00           3.0     500
 30. 24.00           4.0     695
 31. 40.00           5.5     1054
 32. 27.00           3.0     486
 33. 18.00           3.0     442
 34. 62.50           5.5     1249
 35. 53.75           5.0     995
 36. 79.50           5.5     1397
 37. end

.
. /* Each observation is on separate line, with items separated by a space.
>
> We could also read in the data from the Excel file on the class website, using
> this command: */
.
. import excel using "lab 9.xlsx", firstrow clear

.
. /* The firstrow option tells Stata that the firstrow of the Excel file gives the variable names.
>
> As usual, we will add labels to clarify the variable definitions, which is
> more convenient for programming than trying to encode the information in long,
> complicated, hard-to-type variable names. */
.
. label var laborhrs "Labor hours for moving job"

. label var rooms "Number of rooms"

. label var volume "Volume of moving job items (cubic feet)"

.
.
. /* First, we will plot the outcome vs each of the two predictors, adding a
> lowess line to check for a non-linear relationship. */
.
. local ylabel : var label laborhrs

. foreach x in rooms volume {
  2.         local xlabel : var label `x'
  3.         twoway (scatter laborhrs `x') (lowess laborhrs `x'), ///
>                 ytitle("`ylabel'") ylabel(, angle(horizontal)) ///
>                 xtitle("`xlabel'") legend(off) ///
>                 name(lowess_`x', replace)
  4.         }

.
. /* These plots suggest that labor hours may be slightly non-linear in the number of rooms, which could hold if places with more rooms also tend to have bigger rooms and thus more items per room, or require more elaborate packing.  The evidence for non-linearity in the effects of volume seems weaker, primarily driven by the last two observations on the right.
>
> Both relationships are sufficiently linear for the correlations to be meaningful, although they only capture the linear part of the relationship and thus may understate the overall relationship in the presence of non-linearity: */
.
. pwcorr laborhrs rooms volume, sig obs

             | laborhrs    rooms   volume
-------------+---------------------------
    laborhrs |   1.0000
             |
             |       36
             |
       rooms |   0.8638   1.0000
             |   0.0000
             |       36       36
             |
      volume |   0.9430   0.9333   1.0000
             |   0.0000   0.0000
             |       36       36       36
             |

. spearman laborhrs rooms volume, stats(rho obs p)

+-----------------+
|  Key            |
+-----------------+
|   rho           |
|   Number of obs |
|   Sig. level    |
+-----------------+

             | laborhrs    rooms   volume
-------------+---------------------------
    laborhrs |   1.0000
             |       36
             |
             |
       rooms |   0.8884   1.0000
             |       36       36
             |   0.0000
             |
      volume |   0.9292   0.9526   1.0000
             |       36       36       36
             |   0.0000   0.0000
             |

.
. /* Clearly, labor hours is very strongly correlated with both predictors, and
> rank transformation of outcome and predictor has minimal effects on the
> pairwise correlations.  Also the two predictors are highly correlated (i.e.,
> collinear), which would complicate fitting and interpreting models including
> both.
>
> Now Fit a regression predicting labor hours using each of the others one at at
> time. */
.
. foreach x in rooms volume {
  2.         regress laborhrs `x'
  3.         }

      Source |       SS           df       MS      Number of obs   =        36
-------------+----------------------------------   F(1, 34)        =     99.97
       Model |  5799.19478         1  5799.19478   Prob > F        =    0.0000
    Residual |  1972.24272        34  58.0071388   R-squared       =    0.7462
-------------+----------------------------------   Adj R-squared   =    0.7388
       Total |   7771.4375        35  222.041071   Root MSE        =    7.6162

-------------+----------------------------------------------------------------
     laborhrs |      Coef.    Std. Err.       t    P>|t|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
       rooms |    11.13846    1.113992     10.00    0.000      8.874559     13.40237
       _cons |   -10.18099    4.115117     -2.47    0.019     -18.54392    -1.818067
-------------+----------------------------------------------------------------

      Source |       SS           df       MS      Number of obs   =        36
-------------+----------------------------------   F(1, 34)        =    272.99
       Model |  6910.71887         1  6910.71887   Prob > F        =    0.0000
    Residual |  860.718633        34  25.3152539   R-squared       =    0.8892
-------------+----------------------------------   Adj R-squared   =    0.8860
       Total |   7771.4375        35  222.041071   Root MSE        =    5.0314

-------------+----------------------------------------------------------------
     laborhrs |      Coef.    Std. Err.       t    P>|t|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
      volume |    .0500803    .0030311     16.52    0.000      .0439204     .0562402
       _cons |    -2.36966    2.073261     -1.14    0.261     -6.583033     1.843713
-------------+----------------------------------------------------------------

.
. /* Look closely at each model output and consider which model is better.
>
> First, R-squared is high for both models, but considerably higher for the model
> using volume to predict labor hours. Recall that R-squared is the square of the
> Pearson correlations we calculated earlier: */
.
. foreach x in rooms volume {
  2.         qui corr laborhrs `x'
  3.         di as result "`=`r(rho)'^2'"
  4.         }
.746219059183119
.8892458913514151

.
. /* Now describe the model from the equation for labor hours as predicted by the
> volumne of items to be moved.  By typing the regress command without arguments,
> we get back the results of the most recently estimated model: */
.
. regress

      Source |       SS           df       MS      Number of obs   =        36
-------------+----------------------------------   F(1, 34)        =    272.99
       Model |  6910.71887         1  6910.71887   Prob > F        =    0.0000
    Residual |  860.718633        34  25.3152539   R-squared       =    0.8892
-------------+----------------------------------   Adj R-squared   =    0.8860
       Total |   7771.4375        35  222.041071   Root MSE        =    5.0314

-------------+----------------------------------------------------------------
     laborhrs |      Coef.    Std. Err.       t    P>|t|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
      volume |    .0500803    .0030311     16.52    0.000      .0439204     .0562402
       _cons |    -2.36966    2.073261     -1.14    0.261     -6.583033     1.843713
-------------+----------------------------------------------------------------

.
. /* The model shows that expected labor hours increase by about 0.05 hours (95%
> confidence interval 0.043, 0.056, p<0.0005) for each additional cubic foot of
> items to be moved. The standard error of the estimate, which is used in
> calculating the t-statistic, the p-value, and the confidence interval, is
> 0.003.
>
> The overall F-statistic for the model (in this case, with a single predictor,
> equal to the the squuare of the t-statistic for volume) is ~273, and again,
> highly statistically significant (P<0.00005).
>
> Notably, the intercept is not-statistically significant, and thus statistically
> indistinguishable from zero.  This may suggest that the essentially none the
> labor hours are used in doing anything other than moving items; it might be
> positive if labor hours also included ancillary tasks such as getting to and
> from the job.
>
> Finally, we can use the model to predict the labor hours for a job involving
> 1200 cubic feet of items: */
.
. lincom _cons + 1200*volume


 ( 1)  1200*volume + _cons = 0


-------------+----------------------------------------------------------------
     laborhrs |      Coef.    Std. Err.       t    P>|t|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
         (1) |    57.72667    1.932593     29.87    0.000      53.79917     61.65417
-------------+----------------------------------------------------------------

.
. /* So the predicted labor hours for moving 1200 cubic feet of items is 57.7
> hours (95% CI 53.8, 61.7 hours).
>
> Returning to the issue of non-linearity, which was not very evident in the
> Lowess smooth of labor hours on volume: we can check this assumption by adding
> a quadratic term to the model.  */
.
. cap drop volume2

. qui gen volume2 = volume^2

. regress laborhrs volume volume2

      Source |       SS           df       MS      Number of obs   =        36
-------------+----------------------------------   F(2, 33)        =    172.91
       Model |  7094.43363         2  3547.21682   Prob > F        =    0.0000
    Residual |  677.003867        33  20.5152687   R-squared       =    0.9129
-------------+----------------------------------   Adj R-squared   =    0.9076
       Total |   7771.4375        35  222.041071   Root MSE        =    4.5294

-------------+----------------------------------------------------------------
     laborhrs |      Coef.    Std. Err.       t    P>|t|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
      volume |    .0148792    .0120755      1.23    0.227     -.0096885     .0394469
     volume2 |    .0000237    7.91e-06      2.99    0.005      7.58e-06     .0000398
       _cons |    8.578305    4.107046      2.09    0.045      .2224565     16.93415
-------------+----------------------------------------------------------------

.
. /* The quadratic term captures any non-linearity in the relationship of volume
> to labor hours, and is statistically significant (P=0.005).  Also R-squared
> increases a bit, from .8882 to .9129. The coefficients for volume in this
> model are not easy to interpret.  That said, with this slightly better fitting
> model, the constant term is positive and also statistically significant,
> suggesting that there may be some fixed labor costs independent of the volume
> of items.  In addition, we can predict the labor hours for that 1200 cubic foot
> job.  Note that we have to use macro syntax to calculate the square of 1200
> within the lincom command, which won't otherwise calculate it on the fly. */
.
. lincom _cons + 1200*volume + `=1200^2'*volume2


 ( 1)  1200*volume + 1440000*volume2 + _cons = 0


-------------+----------------------------------------------------------------
     laborhrs |      Coef.    Std. Err.       t    P>|t|      [95% Con f. Interval]
-------------+----------------------------------------------------------------
         (1) |    60.51209    1.973104     30.67    0.000      56.49778      64.5264
-------------+----------------------------------------------------------------

.
. /* So this model predicts 60.5 labor hours (95% CI 56.5-64.5) for the 1200
> cubic foot job.
>
> Finally, check the model for outliers. */
.
. predict resid, resid

. graph box resid, name(box_resid, replace)

. qnorm resid, name(qq_resid, replace)

.
. /* No reason for concern here.  In Biostat 208, you will learn to check for
> "influential" points: i.e., observations that change the regression coefficients
> by meaningful amounts, and which usually represent big residuals in combination
> with high-leverage predictor values in the tails of the predictor distribution. */
.
. log close
      name:  <unnamed>
       log:  /Users/eric/Work/MB/ATCR/200/Labs/lab 9.smcl
  log type:  smcl
 closed on:  15 Oct 2018, 08:36:27
-------------