Trulli — https://en.wikipedia.org/wiki/Seattle_box

Overview

In this practical you’ll do basic statistics in R. By the end of this practical you will know how to:

Calculate basic descriptive statistics using mean(), median(), table().
Conduct 1 and 2 sample hypothesis tests with t.test(), cor.test()
Calculate regression analyses with glm() and lm()
Explore statistical objects with names(), summary(), print(), predict()
Use sampling functions such as rnorm() to conduct simulations

Datasets

library(tidyverse)
kc_house <- read_csv("https://raw.githubusercontent.com/therbootcamp/BaselRBootcamp_2018July/master/_sessions/_data//baselrbootcamp_data/kc_house.csv")

File	Rows	Columns	Description
kc_house.csv	21613	21	House sale prices for King County between May 2014 and May 2015.

Packages

Package	Installation
`tidyverse`	`install.packages("tidyverse")`
`lubridate`	`install.packages("lubridate")`
`broom`	`install.packages("broom")`
`rsq`	`install.packages("rsq")`

Glossary

Descriptive Statistics

Function	Description
`table()`	Frequency table
`mean(), median(), mode()`	Measures of central tendency
`sd(), range(), var()`	Measures of variability
`max(), min()`	Extreme values
`summary()`	Several summary statistics

Statistical Tests

Function	Hypothesis Test
`t.test()`	One and two sample t-test
`cor.test()`	Correlation test
`glm()`, `lm()`	Generalized linear model and linear model

Sampling Functions

Function	Description	Additional Help
`sample()`	Draw a random sample of values from a vector	`?sample`
`rnorm()`	Draw random values from a Normal distribution	`?rnorm()`
`runif()`	Draw random values from a Uniform distribution	`?runif()`

Examples

# Examples of hypothesis tests on the diamonds -------------

library(tidyverse)
library(broom)
library(rsq)

# First few rows of the diamonds data

diamonds

# Descriptive statistics ------------------------

mean(diamonds$carat)   # What is the mean carat?
median(diamonds$price)   # What is the median price?
max(diamonds$depth)    # What is the maximum depth?
table(diamonds$color)    # How many observations for each color?


# 2-sample t- test ---------------------------

# Q: Is there a difference in the carats of color = E and color = I diamonds?

htest_B <- t.test(formula = carat ~ color,     # DV ~ IV
                   alternative = "two.sided",  # Two-sided test
                   data = diamonds,         # The data
                   subset = color %in% c("E", "I")) # Compare Diet 1 and Diet 2

htest_B  # Print result


# Correlation test -----------------------

# Q: Is there a correlation between carat and price?

htest_C <- cor.test(formula = ~ carat + price,
                    data = diamonds)

htest_C

# A: Yes. r = 0.92, t(53938) = 551.51, p < .001


# Regression ----------------------------


# Q: Create regression equation predicting price by carat, depth, table, and x

price_glm <- glm(formula = price ~ carat + depth + table + x,
                  data = diamonds)

# Print coefficients
price_glm$coefficients

# Tidy version
tidy(price_glm)

# Extract R-Squared
rsq(price_glm)

# -----
# Simulation
# ------

# 100 random samples from a normal distribution with mean = 0, sd = 1
samp_A <- rnorm(n = 100, mean = 0, sd = 1)

# 100 random samples from a Uniform distribution with bounds at 0, 10
samp_B <- runif(n = 100, min = 0, max = 10)

# Calculate descriptives
mean(samp_A)
sd(samp_A)

mean(samp_B)
sd(samp_B)

# Combine samples (plus tw new ones) in a tibble

my_sim <- tibble(A = samp_A,
                 B = samp_B,
                 C = rnorm(n = 100, mean = 0, sd = 1),
                 error = rnorm(n = 100, mean = 5, sd = 10))

# Add y, a linear function of A and B to my_sim
my_sim <- my_sim %>%
  mutate(y = 3 * A -8 * B + error)

# Regress y on A, B and C
my_glm <- glm(y ~ A + B + C,
              data = my_sim)

# Look at results!
tidy(my_glm)

Tasks

A - Getting setup

Open your baselrbootcamp R project. It should already have the folders 1_Data and 2_Code. Make sure that all of the data files listed above are contained in the 1_Data folder

# Done!

Open a new R script and save it as a new file called statistics_practical.R in the 2_Code folder. At the top of the script, using comments, write your name and the date. The, load the set of packages listed above with library().
For this practical, we’ll use the kc_house.csv data. This dataset contains house sale prices for King County, Washington. It includes homes sold between May 2014 and May 2015. Using the following template, load the data into R and store it as a new object called kc_house.

kc_house <- read_csv(file = "XX")

kc_house <- read_csv(file = "1_Data/kc_house.csv")

Using print(), summary(), and head(), explore the data to make sure it was loaded correctly.

kc_house

# A tibble: 21,613 x 21
   id         date                  price bedrooms bathrooms sqft_living
   <chr>      <dttm>                <dbl>    <int>     <dbl>       <int>
 1 7129300520 2014-10-13 00:00:00  221900        3      1           1180
 2 6414100192 2014-12-09 00:00:00  538000        3      2.25        2570
 3 5631500400 2015-02-25 00:00:00  180000        2      1            770
 4 2487200875 2014-12-09 00:00:00  604000        4      3           1960
 5 1954400510 2015-02-18 00:00:00  510000        3      2           1680
 6 7237550310 2014-05-12 00:00:00 1225000        4      4.5         5420
 7 1321400060 2014-06-27 00:00:00  257500        3      2.25        1715
 8 2008000270 2015-01-15 00:00:00  291850        3      1.5         1060
 9 2414600126 2015-04-15 00:00:00  229500        3      1           1780
10 3793500160 2015-03-12 00:00:00  323000        3      2.5         1890
# ... with 21,603 more rows, and 15 more variables: sqft_lot <int>,
#   floors <dbl>, waterfront <int>, view <int>, condition <int>,
#   grade <int>, sqft_above <int>, sqft_basement <int>, yr_built <int>,
#   yr_renovated <int>, zipcode <int>, lat <dbl>, long <dbl>,
#   sqft_living15 <int>, sqft_lot15 <int>

summary(kc_house)

      id                 date                         price        
 Length:21613       Min.   :2014-05-02 00:00:00   Min.   :  75000  
 Class :character   1st Qu.:2014-07-22 00:00:00   1st Qu.: 321950  
 Mode  :character   Median :2014-10-16 00:00:00   Median : 450000  
                    Mean   :2014-10-29 04:38:01   Mean   : 540088  
                    3rd Qu.:2015-02-17 00:00:00   3rd Qu.: 645000  
                    Max.   :2015-05-27 00:00:00   Max.   :7700000  
    bedrooms      bathrooms     sqft_living       sqft_lot      
 Min.   : 0.0   Min.   :0.00   Min.   :  290   Min.   :    520  
 1st Qu.: 3.0   1st Qu.:1.75   1st Qu.: 1427   1st Qu.:   5040  
 Median : 3.0   Median :2.25   Median : 1910   Median :   7618  
 Mean   : 3.4   Mean   :2.11   Mean   : 2080   Mean   :  15107  
 3rd Qu.: 4.0   3rd Qu.:2.50   3rd Qu.: 2550   3rd Qu.:  10688  
 Max.   :33.0   Max.   :8.00   Max.   :13540   Max.   :1651359  
     floors       waterfront         view        condition   
 Min.   :1.00   Min.   :0.000   Min.   :0.00   Min.   :1.00  
 1st Qu.:1.00   1st Qu.:0.000   1st Qu.:0.00   1st Qu.:3.00  
 Median :1.50   Median :0.000   Median :0.00   Median :3.00  
 Mean   :1.49   Mean   :0.008   Mean   :0.23   Mean   :3.41  
 3rd Qu.:2.00   3rd Qu.:0.000   3rd Qu.:0.00   3rd Qu.:4.00  
 Max.   :3.50   Max.   :1.000   Max.   :4.00   Max.   :5.00  
     grade         sqft_above   sqft_basement     yr_built   
 Min.   : 1.00   Min.   : 290   Min.   :   0   Min.   :1900  
 1st Qu.: 7.00   1st Qu.:1190   1st Qu.:   0   1st Qu.:1951  
 Median : 7.00   Median :1560   Median :   0   Median :1975  
 Mean   : 7.66   Mean   :1788   Mean   : 292   Mean   :1971  
 3rd Qu.: 8.00   3rd Qu.:2210   3rd Qu.: 560   3rd Qu.:1997  
 Max.   :13.00   Max.   :9410   Max.   :4820   Max.   :2015  
  yr_renovated     zipcode           lat            long     
 Min.   :   0   Min.   :98001   Min.   :47.2   Min.   :-122  
 1st Qu.:   0   1st Qu.:98033   1st Qu.:47.5   1st Qu.:-122  
 Median :   0   Median :98065   Median :47.6   Median :-122  
 Mean   :  84   Mean   :98078   Mean   :47.6   Mean   :-122  
 3rd Qu.:   0   3rd Qu.:98118   3rd Qu.:47.7   3rd Qu.:-122  
 Max.   :2015   Max.   :98199   Max.   :47.8   Max.   :-121  
 sqft_living15    sqft_lot15    
 Min.   : 399   Min.   :   651  
 1st Qu.:1490   1st Qu.:  5100  
 Median :1840   Median :  7620  
 Mean   :1987   Mean   : 12768  
 3rd Qu.:2360   3rd Qu.: 10083  
 Max.   :6210   Max.   :871200

head(kc_house)

# A tibble: 6 x 21
  id    date                 price bedrooms bathrooms sqft_living sqft_lot
  <chr> <dttm>               <dbl>    <int>     <dbl>       <int>    <int>
1 7129… 2014-10-13 00:00:00 2.22e5        3      1           1180     5650
2 6414… 2014-12-09 00:00:00 5.38e5        3      2.25        2570     7242
3 5631… 2015-02-25 00:00:00 1.80e5        2      1            770    10000
4 2487… 2014-12-09 00:00:00 6.04e5        4      3           1960     5000
5 1954… 2015-02-18 00:00:00 5.10e5        3      2           1680     8080
6 7237… 2014-05-12 00:00:00 1.23e6        4      4.5         5420   101930
# ... with 14 more variables: floors <dbl>, waterfront <int>, view <int>,
#   condition <int>, grade <int>, sqft_above <int>, sqft_basement <int>,
#   yr_built <int>, yr_renovated <int>, zipcode <int>, lat <dbl>,
#   long <dbl>, sqft_living15 <int>, sqft_lot15 <int>

B - Descriptive statistics

What was the mean price (price) of all houses?

mean(kc_house$price)

[1] 540088

What was the median price (price) of all houses?

median(kc_house$price)

[1] 450000

What was the standard deviation price (price) of all houses?

sd(kc_house$price)

[1] 367127

What were the minimum and maximum prices? Try using min(), max(), and range()

range(kc_house$price)

[1]   75000 7700000

How many houses were on the waterfront (waterfront) and how many were not? (Hint: use table())

table(kc_house$waterfront)


    0     1 
21450   163

C - T tests with t.test()

Do houses on the waterfront sell for more than those not on the waterfront? Answer this by creating an object called water_htest which contains the result of the appropriate t-test using t.test(). Use the following template

water_htest <- t.test(formula = XX ~ XX,
                      data = XX)

water_htest <- t.test(formula = price ~ waterfront,
                      data = kc_house)

Print your water_htest object to see summary results

water_htest


    Welch Two Sample t-test

data:  price by waterfront
t = -10, df = 200, p-value <2e-16
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -1303662  -956963
sample estimates:
mean in group 0 mean in group 1 
         531564         1661876

Apply the summary() function to water_htest. Do you see anything new or important?

summary(water_htest)

            Length Class  Mode     
statistic   1      -none- numeric  
parameter   1      -none- numeric  
p.value     1      -none- numeric  
conf.int    2      -none- numeric  
estimate    2      -none- numeric  
null.value  1      -none- numeric  
alternative 1      -none- character
method      1      -none- character
data.name   1      -none- character

Look at the names of the water_htest with names(). What are the named elements of the object?

names(water_htest)

[1] "statistic"   "parameter"   "p.value"     "conf.int"    "estimate"   
[6] "null.value"  "alternative" "method"      "data.name"

Using the $ operator, print the exact test statistic of the t-test.

water_htest$statistic

    t 
-12.9

Again using the $ operator, print the exact p-value of the test.

water_htest$p.value

[1] 1.38e-26

Use the tidy() function to return a dataframe containing the main results of the test.

tidy(water_htest)

  estimate estimate1 estimate2 statistic  p.value parameter conf.low
1 -1130312    531564   1661876     -12.9 1.38e-26       162 -1303662
  conf.high                  method alternative
1   -956963 Welch Two Sample t-test   two.sided

Do houses on the waterfront also tend to have more (or less) living space (sqft_living)? Answer this by conducting the appropriate t-test and going through the same steps you did in the previous question

D - Correlation tests with cor.test()

Do houses built in later years tend to sell for more than older houses? Answer this by creating an object called time_htest which contains the result of the appropriate correlation test using cor.test() and the columns yr_built and price. Use the following template:

# Note: cor.terst() is the only function (we know of)
#  that uses formulas in the strange format
#   formula = ~ X + Y instead of formula = Y ~ X

time_htest <- cor.test(formula = ~ XX + XX,
                       data = XX)

time_htest <- cor.test(formula = ~ yr_built + price,
                       data = kc_house)

Print your time_htest object to see summary results

time_htest


    Pearson's product-moment correlation

data:  yr_built and price
t = 8, df = 20000, p-value = 2e-15
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.0407 0.0673
sample estimates:
  cor 
0.054

Apply the summary() function to time_htest. Do you see anything new or important?

summary(time_htest)

            Length Class  Mode     
statistic   1      -none- numeric  
parameter   1      -none- numeric  
p.value     1      -none- numeric  
estimate    1      -none- numeric  
null.value  1      -none- numeric  
alternative 1      -none- character
method      1      -none- character
data.name   1      -none- character
conf.int    2      -none- numeric

Look at the names of the time_htest with names(). What are the named elements of the object?

names(time_htest)

[1] "statistic"   "parameter"   "p.value"     "estimate"    "null.value" 
[6] "alternative" "method"      "data.name"   "conf.int"

Using the $ operator, print the exact correlation of the test (Hint: it’s called “estimate”)

time_htest$estimate

  cor 
0.054

Again using the $ operator, print the exact p-value of the test.

time_htest$p.value

[1] 1.93e-15

Use the tidy() function to return a dataframe containing the main results of the test.

tidy(time_htest)

  estimate statistic  p.value parameter conf.low conf.high
1    0.054      7.95 1.93e-15     21611   0.0407    0.0673
                                method alternative
1 Pearson's product-moment correlation   two.sided

Do houses in in better overall condition condition tend to sell at higher prices than those in worse condition? Answer this by conducting the appropriate correlation test and going through the same steps you did in the previous question.

condition_htest <- cor.test(formula = ~ condition + price,
                       data = kc_house)

condition_htest


    Pearson's product-moment correlation

data:  condition and price
t = 5, df = 20000, p-value = 9e-08
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.0230 0.0497
sample estimates:
   cor 
0.0364

E - Linear Regression with lm() and glm()

In a previous question, you used a correlation test to see if houses built in later years yr_built tend to sell for more than older houses. Now let’s repeat the same analysis using a regression analysis with lm(). Do to this, use the following template:

time_lm <- lm(formula = XX ~ XX,
                 data = XX)

time_lm <- lm(formula = price ~ yr_built,
                 data = kc_house)

Print your time_lm object to see the main results.

time_lm


Call:
lm(formula = price ~ yr_built, data = kc_house)

Coefficients:
(Intercept)     yr_built  
    -790478          675

Use the tidy() function (from the broom package) to return a dataframe containing the main results of the test.

tidy(time_lm)

         term estimate std.error statistic  p.value
1 (Intercept)  -790478  167350.2     -4.72 2.33e-06
2    yr_built      675      84.9      7.95 1.93e-15

How does the p-value of this test compare to what you got in the previous question when you used cor.test()?

# Same!

Which of the variables bedrooms, bathrooms, and sqft_living predict housing price? Answer this by conducting the appropriate regression analysis using lm() and assigning the result to an object price_lm. Use the following template:

price_lm <- lm(formula = XX ~ XX + XX + XX,
               data = XX)

price_lm <- lm(formula = price ~ bedrooms + bathrooms + sqft_living,
               data = kc_house)

Print your price_lm object to see the main results.

price_lm


Call:
lm(formula = price ~ bedrooms + bathrooms + sqft_living, data = kc_house)

Coefficients:
(Intercept)     bedrooms    bathrooms  sqft_living  
      74847       -57861         7933          309

Apply the summary() function to price_lm. Do you see anything new or important?

summary(price_lm)


Call:
lm(formula = price ~ bedrooms + bathrooms + sqft_living, data = kc_house)

Residuals:
     Min       1Q   Median       3Q      Max 
-1642456  -144268   -22821   102461  4180678 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)  74847.14    6913.67   10.83   <2e-16 ***
bedrooms    -57860.89    2334.61  -24.78   <2e-16 ***
bathrooms     7932.71    3510.56    2.26    0.024 *  
sqft_living    309.39       3.09  100.23   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 258000 on 21609 degrees of freedom
Multiple R-squared:  0.507, Adjusted R-squared:  0.507 
F-statistic: 7.41e+03 on 3 and 21609 DF,  p-value: <2e-16

Look at the names of the price_lm with names(). What are the named elements of the object?

names(price_lm)

 [1] "coefficients"  "residuals"     "effects"       "rank"         
 [5] "fitted.values" "assign"        "qr"            "df.residual"  
 [9] "xlevels"       "call"          "terms"         "model"

Using the $ operator, print a vector of the estimated coefficients of the analysis.

price_lm$coefficients

(Intercept)    bedrooms   bathrooms sqft_living 
      74847      -57861        7933         309

Use the tidy() function (from the broom package) to return a dataframe containing the main results of the test.

tidy(price_lm)

         term estimate std.error statistic   p.value
1 (Intercept)    74847   6913.67     10.83  3.05e-27
2    bedrooms   -57861   2334.61    -24.78 9.81e-134
3   bathrooms     7933   3510.56      2.26  2.39e-02
4 sqft_living      309      3.09    100.23  0.00e+00

Use the rsq() function (from the rsq package) to obtain the R-squared value. Does this match what you saw in your previous outputs?

rsq(price_lm)

[1] 0.507

Using the following template, create a new model called everything_lm that predicts housing prices based on all predictors in kc_house except for id and date (id is meaningless and date shouldn’t matter). Use the following template. Note that to include all predictors, we’ll use the formula = y ~. shortcut. We’ll also remove id and date from the dataset using select() before running the analysis.

everything_lm <- lm(formula = XX ~.,
                    data = kc_house %>% 
                             select(-id, -date))

everything_lm <- lm(formula = price ~.,
                    data = kc_house %>% 
                             select(-id, -date))

Print your everything_lm object to see the main results.

everything_lm


Call:
lm(formula = price ~ ., data = kc_house %>% select(-id, -date))

Coefficients:
  (Intercept)       bedrooms      bathrooms    sqft_living       sqft_lot  
     6.69e+06      -3.58e+04       4.11e+04       1.50e+02       1.29e-01  
       floors     waterfront           view      condition          grade  
     6.69e+03       5.83e+05       5.29e+04       2.64e+04       9.59e+04  
   sqft_above  sqft_basement       yr_built   yr_renovated        zipcode  
     3.11e+01             NA      -2.62e+03       1.98e+01      -5.82e+02  
          lat           long  sqft_living15     sqft_lot15  
     6.03e+05      -2.15e+05       2.17e+01      -3.83e-01

Apply the summary() function to everything_lm. Do you see anything new or important?

summary(everything_lm)


Call:
lm(formula = price ~ ., data = kc_house %>% select(-id, -date))

Residuals:
     Min       1Q   Median       3Q      Max 
-1291725   -99229    -9739    77583  4333222 

Coefficients: (1 not defined because of singularities)
               Estimate Std. Error t value Pr(>|t|)    
(Intercept)    6.69e+06   2.93e+06    2.28   0.0225 *  
bedrooms      -3.58e+04   1.89e+03  -18.91  < 2e-16 ***
bathrooms      4.11e+04   3.25e+03   12.65  < 2e-16 ***
sqft_living    1.50e+02   4.39e+00   34.23  < 2e-16 ***
sqft_lot       1.29e-01   4.79e-02    2.68   0.0073 ** 
floors         6.69e+03   3.60e+03    1.86   0.0628 .  
waterfront     5.83e+05   1.74e+04   33.58  < 2e-16 ***
view           5.29e+04   2.14e+03   24.71  < 2e-16 ***
condition      2.64e+04   2.35e+03   11.22  < 2e-16 ***
grade          9.59e+04   2.15e+03   44.54  < 2e-16 ***
sqft_above     3.11e+01   4.36e+00    7.14  9.7e-13 ***
sqft_basement        NA         NA      NA       NA    
yr_built      -2.62e+03   7.27e+01  -36.06  < 2e-16 ***
yr_renovated   1.98e+01   3.66e+00    5.42  6.0e-08 ***
zipcode       -5.82e+02   3.30e+01  -17.66  < 2e-16 ***
lat            6.03e+05   1.07e+04   56.15  < 2e-16 ***
long          -2.15e+05   1.31e+04  -16.35  < 2e-16 ***
sqft_living15  2.17e+01   3.45e+00    6.29  3.3e-10 ***
sqft_lot15    -3.83e-01   7.33e-02   -5.22  1.8e-07 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 201000 on 21595 degrees of freedom
Multiple R-squared:   0.7,  Adjusted R-squared:   0.7 
F-statistic: 2.96e+03 on 17 and 21595 DF,  p-value: <2e-16

Using the $ operator, print a vector of the estimated coefficients of the analysis. Are the beta values for bedrooms, bathrooms, and sqft_living different from what you got in your previous analysis price_lm?

everything_lm$coefficients

  (Intercept)      bedrooms     bathrooms   sqft_living      sqft_lot 
     6.69e+06     -3.58e+04      4.11e+04      1.50e+02      1.29e-01 
       floors    waterfront          view     condition         grade 
     6.69e+03      5.83e+05      5.29e+04      2.64e+04      9.59e+04 
   sqft_above sqft_basement      yr_built  yr_renovated       zipcode 
     3.11e+01            NA     -2.62e+03      1.98e+01     -5.82e+02 
          lat          long sqft_living15    sqft_lot15 
     6.03e+05     -2.15e+05      2.17e+01     -3.83e-01

Use the tidy() function to return a dataframe containing the main results of the test.

tidy(everything_lm)

            term  estimate std.error statistic   p.value
1    (Intercept)  6.69e+06  2.93e+06      2.28  2.25e-02
2       bedrooms -3.58e+04  1.89e+03    -18.91  4.46e-79
3      bathrooms  4.11e+04  3.25e+03     12.65  1.60e-36
4    sqft_living  1.50e+02  4.39e+00     34.23 4.47e-250
5       sqft_lot  1.29e-01  4.79e-02      2.68  7.29e-03
6         floors  6.69e+03  3.60e+03      1.86  6.28e-02
7     waterfront  5.83e+05  1.74e+04     33.58 5.01e-241
8           view  5.29e+04  2.14e+03     24.71 6.54e-133
9      condition  2.64e+04  2.35e+03     11.22  3.87e-29
10         grade  9.59e+04  2.15e+03     44.54  0.00e+00
11    sqft_above  3.11e+01  4.36e+00      7.14  9.71e-13
12      yr_built -2.62e+03  7.27e+01    -36.06 1.39e-276
13  yr_renovated  1.98e+01  3.66e+00      5.42  6.03e-08
14       zipcode -5.82e+02  3.30e+01    -17.66  2.78e-69
15           lat  6.03e+05  1.07e+04     56.15  0.00e+00
16          long -2.15e+05  1.31e+04    -16.35  1.01e-59
17 sqft_living15  2.17e+01  3.45e+00      6.29  3.26e-10
18    sqft_lot15 -3.83e-01  7.33e-02     -5.22  1.78e-07

Use the rsq() function (from the rsq package) to obtain the R-squared value. How does this R-squared value compare to what you got in your previous regression price_lm?

rsq(everything_lm)

[1] 0.7

How well did the everything_lm model fit the actual housing prices? We can answer this by calculating the average difference between the fitted values and the true values directly. Using the following template, make this calculation. What do you find is the mean absolute difference between the fitted housing prices and the true prices?

# True housing prices
prices_true <- kc_house$XX

# Model fits (fitted.values)
prices_fitted <- everything_lm$XX

# Calculate absolute error between fitted and true values
abs_error <- abs(XX - XX)

# Calculate mean absolute error
mae <- mean(XX)

# Print the result!
mae

# True housing prices
prices_true <- kc_house$price

# Model fits (fitted.values)
prices_fitted <- everything_lm$fitted.values

# Calculate absolute error between fitted and true values
abs_error <- abs(prices_true - prices_fitted)

# Calculate mean absolute error
mae <- mean(abs_error)

# Print the result!
mae

[1] 125923

Using the following template, create a scatterplot showing the relationship between the fitted values of the everything_lm object and the true prices.

# Create dataframe containing fitted and true prices
prices_results <- tibble(truth = XX,
                         fitted = XX)

# Create scatterplot
ggplot(data = prices_results,
       aes(x = fitted, y = truth)) + 
  geom_point(alpha = .1) +
  geom_smooth(method = 'lm') +
  labs(title = "Fitted versus Predicted housing prices",
       subtitle = paste0("Mean Absolute Error = ", mae),
       caption = "Source: King County housing price Kaggle dataset",
       x = "Fitted housing price values",
       y = 'True values') + 
  theme_minimal()

# Create dataframe containing fitted and true prices
prices_results <- tibble(truth = kc_house$price,
                         fitted = everything_lm$fitted.values)

# Create scatterplot
ggplot(data = prices_results,
       aes(x = fitted, y = truth)) + 
  geom_point(alpha = .1) +
  geom_smooth(method = 'lm') +
  labs(title = "Fitted versus Predicted housing prices",
       subtitle = paste0("Mean Absolute Error = ", mae),
       caption = "Source: King County housing price Kaggle dataset",
       x = "Fitted housing price values",
       y = 'True values') + 
  theme_minimal()

Statistics

Basel R Bootcamp July 2018www.therbootcamp.com@therbootcamp

Overview

Datasets

Packages

Glossary

Examples

Tasks

A - Getting setup

B - Descriptive statistics

C - T tests with t.test()

D - Correlation tests with cor.test()

E - Linear Regression with lm() and glm()

Basel R Bootcamp July 2018
www.therbootcamp.com
@therbootcamp