Data Analysis
TTTT
Sep 14, 2020
1. Read Data
t = "F:\\TAI LIEU THAM KHAO\\Thesis Thao\\du lieu.csv"
DL = read.csv(t, header = T)
head(DL)## Temp.pav Hour Depth Method Temp.env
## 1 31.00 24 4 1 32.61
## 2 34.75 24 4 2 32.61
## 3 30.25 24 4 3 32.61
## 4 30.75 24 8 1 32.61
## 5 35.75 24 8 2 32.61
## 6 31.50 24 8 3 32.61fix(DL)
names(DL)## [1] "Temp.pav" "Hour" "Depth" "Method" "Temp.env"str(DL)## 'data.frame': 3159 obs. of 5 variables:
## $ Temp.pav: num 31 34.8 30.2 30.8 35.8 ...
## $ Hour : num 24 24 24 24 24 24 24 24 24 1 ...
## $ Depth : num 4 4 4 8 8 8 12 12 12 4 ...
## $ Method : num 1 2 3 1 2 3 1 2 3 1 ...
## $ Temp.env: Factor w/ 103 levels " - ","26.9",..: 56 56 56 56 56 56 56 56 56 18 ...dim(DL)## [1] 3159 5summary(DL)## Temp.pav Hour Depth Method Temp.env
## Min. :29.25 Min. : 1.0 Min. : 4 Min. :1 29.2 : 90
## 1st Qu.:34.25 1st Qu.: 6.0 1st Qu.: 4 1st Qu.:1 30.1 : 81
## Median :36.75 Median :13.0 Median : 8 Median :2 30.9 : 81
## Mean :37.63 Mean :12.5 Mean : 8 Mean :2 28.3 : 72
## 3rd Qu.:40.00 3rd Qu.:19.0 3rd Qu.:12 3rd Qu.:3 28.8 : 72
## Max. :57.00 Max. :24.0 Max. :12 Max. :3 31 : 72
## (Other):2691attach(DL)
DL = within(DL, {
Temp.pav = as.numeric(Temp.pav)
Hour = as.numeric(Hour)
Depth = as.numeric(Depth)
Method = factor(Method, labels = c("VD", "NT", "BB"))
Temp.env = as.numeric(Temp.env)
} )
head(DL)## Temp.pav Hour Depth Method Temp.env
## 1 31.00 24 4 VD 56
## 2 34.75 24 4 NT 56
## 3 30.25 24 4 BB 56
## 4 30.75 24 8 VD 56
## 5 35.75 24 8 NT 56
## 6 31.50 24 8 BB 56str(DL)## 'data.frame': 3159 obs. of 5 variables:
## $ Temp.pav: num 31 34.8 30.2 30.8 35.8 ...
## $ Hour : num 24 24 24 24 24 24 24 24 24 1 ...
## $ Depth : num 4 4 4 8 8 8 12 12 12 4 ...
## $ Method : Factor w/ 3 levels "VD","NT","BB": 1 2 3 1 2 3 1 2 3 1 ...
## $ Temp.env: num 56 56 56 56 56 56 56 56 56 18 ...2. Descriptive analyse
# Descriptive analyse with package "table 1" - sort by Bus variable
## install.packages("table1", dependencies = T)
library(table1)##
## Attaching package: 'table1'## The following objects are masked from 'package:base':
##
## units, units<-table1 (~ Temp.pav + Hour + Depth + Temp.env| Method , data = DL)| VD (n=1053) |
NT (n=1053) |
BB (n=1053) |
Overall (n=3159) |
|
|---|---|---|---|---|
| Temp.pav | ||||
| Mean (SD) | 37.0 (4.40) | 40.6 (5.31) | 35.3 (2.86) | 37.6 (4.84) |
| Median [Min, Max] | 36.3 [29.3, 51.5] | 39.3 [31.5, 57.0] | 35.3 [29.3, 42.8] | 36.8 [29.3, 57.0] |
| Hour | ||||
| Mean (SD) | 12.5 (7.00) | 12.5 (7.00) | 12.5 (7.00) | 12.5 (7.00) |
| Median [Min, Max] | 13.0 [1.00, 24.0] | 13.0 [1.00, 24.0] | 13.0 [1.00, 24.0] | 13.0 [1.00, 24.0] |
| Depth | ||||
| Mean (SD) | 8.00 (3.27) | 8.00 (3.27) | 8.00 (3.27) | 8.00 (3.27) |
| Median [Min, Max] | 8.00 [4.00, 12.0] | 8.00 [4.00, 12.0] | 8.00 [4.00, 12.0] | 8.00 [4.00, 12.0] |
| Temp.env | ||||
| Mean (SD) | 46.7 (25.8) | 46.7 (25.8) | 46.7 (25.8) | 46.7 (25.8) |
| Median [Min, Max] | 42.0 [1.00, 103] | 42.0 [1.00, 103] | 42.0 [1.00, 103] | 42.0 [1.00, 103] |
table1 (~ Temp.pav + Hour + Depth + Temp.env + Method , data = DL)| Overall (n=3159) |
|
|---|---|
| Temp.pav | |
| Mean (SD) | 37.6 (4.84) |
| Median [Min, Max] | 36.8 [29.3, 57.0] |
| Hour | |
| Mean (SD) | 12.5 (7.00) |
| Median [Min, Max] | 13.0 [1.00, 24.0] |
| Depth | |
| Mean (SD) | 8.00 (3.27) |
| Median [Min, Max] | 8.00 [4.00, 12.0] |
| Temp.env | |
| Mean (SD) | 46.7 (25.8) |
| Median [Min, Max] | 42.0 [1.00, 103] |
| Method | |
| VD | 1053 (33.3%) |
| NT | 1053 (33.3%) |
| BB | 1053 (33.3%) |
# Descriptive analyse with package "compareGroups" - sort by Bus variable (hon table1 la cung cap them tri so p. Doi voi bien phan loai, tri so p cung duoc tinh theo Chi-Square)
library(compareGroups)## Loading required package: SNPassoc## Loading required package: haplo.stats## Loading required package: survival## Loading required package: mvtnorm## Loading required package: parallel## Registered S3 method overwritten by 'SNPassoc':
## method from
## summary.haplo.glm haplo.statst = compareGroups(Method ~ Temp.pav + Hour + Depth + Temp.env, data = DL)
createTable(t)##
## --------Summary descriptives table by 'Method'---------
##
## ______________________________________________________
## VD NT BB p.overall
## N=1053 N=1053 N=1053
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
## Temp.pav 37.0 (4.40) 40.6 (5.31) 35.3 (2.86) <0.001
## Hour 12.5 (7.00) 12.5 (7.00) 12.5 (7.00) 1.000
## Depth 8.00 (3.27) 8.00 (3.27) 8.00 (3.27) 1.000
## Temp.env 46.7 (25.8) 46.7 (25.8) 46.7 (25.8) 1.000
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯3. Descriptive analyse by graphs
library(magrittr)
library(tidyverse)## -- Attaching packages -------------------------------------------------------------- tidyverse 1.3.0 --## v ggplot2 3.2.1 v purrr 0.3.3
## v tibble 2.1.3 v dplyr 0.8.3
## v tidyr 1.0.0 v stringr 1.4.0
## v readr 1.3.1 v forcats 0.4.0## -- Conflicts ----------------------------------------------------------------- tidyverse_conflicts() --
## x tidyr::extract() masks magrittr::extract()
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
## x purrr::set_names() masks magrittr::set_names()library(ggplot2)
library(car)## Loading required package: carData## Registered S3 methods overwritten by 'car':
## method from
## influence.merMod lme4
## cooks.distance.influence.merMod lme4
## dfbeta.influence.merMod lme4
## dfbetas.influence.merMod lme4##
## Attaching package: 'car'## The following object is masked from 'package:dplyr':
##
## recode## The following object is masked from 'package:purrr':
##
## some# Correlation
## All variables
library(psych)##
## Attaching package: 'psych'## The following object is masked from 'package:car':
##
## logit## The following objects are masked from 'package:ggplot2':
##
## %+%, alphaattach(DL)## The following objects are masked from DL (pos = 23):
##
## Depth, Hour, Method, Temp.env, Temp.pavvars = cbind(Temp.pav, Hour, Depth, Method, Temp.env)
pairs.panels(vars)
pairs.panels(DL)
## Only continuous variables
vars1 = cbind(Temp.pav, Hour, Depth, Temp.env)
pairs.panels(vars1)
# Statistic Analysis by boxplot (continuous variales)
DL %>%
group_by(Temp.pav, Method) %>%
count() %>%
ggplot(aes(x = Method, y = Temp.pav, fill = Method)) +
geom_boxplot() +
xlab("Method") +
ylab("Temperature of pavement") +
ggtitle("Boxplot of Temperature of Pavement ~ Method")
DL %>%
group_by(Temp.env, Method) %>%
count() %>%
ggplot(aes(x = Method, y = Temp.env, fill = Method)) +
geom_boxplot() +
xlab("Method") +
ylab("Temperature of environment") +
ggtitle("Boxplot of Temperature of environment ~ Method")
DL %>%
group_by(Hour, Method) %>%
count() %>%
ggplot(aes(x = Method, y = Hour, fill = Method)) +
geom_boxplot() +
xlab("Method") +
ylab("Hour of day") +
ggtitle("Boxplot of Hour of date ~ Method")
DL %>%
group_by(Depth, Method) %>%
count() %>%
ggplot(aes(x = Method, y = Depth, fill = Method)) +
geom_boxplot() +
xlab("Method") +
ylab("Depth") +
ggtitle("Boxplot of Depth ~ Method")
# Statistic Analysis by Histogram (Countinuous Variables)
ggplot(DL, aes(x=Temp.pav, fill = Method, color = Method)) +
geom_histogram (position = "dodge") +
xlab("Temperature of pavement") +
ylab("Count") +
ggtitle("Histogram of Temperature of pavement")## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
ggplot(DL, aes(x=Temp.env, fill = Method, color = Method)) +
geom_histogram (position = "dodge") +
xlab("Temperature of environment") +
ylab("Count") +
ggtitle("Histogram of Temperature of environment")## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
# Statistic Analysis by geom_point - Correlation with a couple of variables (Tuong quan don gian = Tan xa)
ggplot(DL, aes(x = Temp.env, y = Temp.pav, fill = Method)) +
geom_point(aes(col = Method, size = 0.000001)) +
geom_smooth(method = "lm", formula = y ~ x + I(x^2)) +
xlab("Temperature of environment") +
ylab("Temperature of pavement") +
ggtitle("Temperature of Pavement ~ Temperature of environment")
ggplot(DL, aes(x = Temp.env, y = Hour, fill = Method)) +
geom_point(aes(col = Method, size = 0.000001)) +
geom_smooth(method = "lm") +
xlab("Temperature of environment") +
ylab("Hour") +
ggtitle("Hour ~ Temperature of environment")
ggplot(DL, aes(x = Hour, y = Temp.pav)) +
geom_point(aes(col = Method, size = 0.000001)) +
xlab("Hour") +
ylab("Temperature of pavement") +
ggtitle("Temperature of Pavement ~ Hour")
ggplot(DL, aes(x = Depth, y = Temp.pav)) +
geom_point(aes(col = Method, size = 0.000001)) +
xlab("Depth") +
ylab("Temperature of pavement") +
ggtitle("Temperature of Pavement ~ Depth")
library(car)
attach(DL)## The following objects are masked from DL (pos = 3):
##
## Depth, Hour, Method, Temp.env, Temp.pav## The following objects are masked from DL (pos = 24):
##
## Depth, Hour, Method, Temp.env, Temp.pavscatterplot(Temp.pav ~ Temp.env)
scatterplot(Temp.pav ~ Temp.env, main = "Temperature of pavement ~ Temperature of environment", xlab = "Temperature of environment", ylab = "Temperature of pavement")
Test the differences of Temperature of pavement between different methoths
# F-test (T-test voi 2 nhom)
## Dung anova
av <- aov(Temp.pav ~ as.factor(Method))
av## Call:
## aov(formula = Temp.pav ~ as.factor(Method))
##
## Terms:
## as.factor(Method) Residuals
## Sum of Squares 15288.70 58643.85
## Deg. of Freedom 2 3156
##
## Residual standard error: 4.31065
## Estimated effects may be unbalancedsummary(av)## Df Sum Sq Mean Sq F value Pr(>F)
## as.factor(Method) 2 15289 7644 411.4 <2e-16 ***
## Residuals 3156 58644 19
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1TukeyHSD(av) # Dieu chinh p theo pp Tukey## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Temp.pav ~ as.factor(Method))
##
## $`as.factor(Method)`
## diff lwr upr p adj
## NT-VD 3.520418 3.079913 3.960923 0
## BB-VD -1.773029 -2.213535 -1.332524 0
## BB-NT -5.293447 -5.733952 -4.852942 0## Dung pairwise.t.test va dieu chinh p d ss nhieu nhom
pairwise.t.test(Temp.pav, Method, p.adjust.method = "bonferroni")##
## Pairwise comparisons using t tests with pooled SD
##
## data: Temp.pav and Method
##
## VD NT
## NT <2e-16 -
## BB <2e-16 <2e-16
##
## P value adjustment method: bonferroniRegression Model
# Simple linear regression models
summary(lm(Temp.pav ~ Hour, data = DL))##
## Call:
## lm(formula = Temp.pav ~ Hour, data = DL)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.6184 -3.2289 -0.7289 2.0503 19.2742
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 35.1936 0.1691 208.14 <2e-16 ***
## Hour 0.1948 0.0118 16.51 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.643 on 3157 degrees of freedom
## Multiple R-squared: 0.07947, Adjusted R-squared: 0.07918
## F-statistic: 272.5 on 1 and 3157 DF, p-value: < 2.2e-16summary(lm(Temp.pav ~ Depth, data = DL))##
## Call:
## lm(formula = Temp.pav ~ Depth, data = DL)
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.2452 -3.3792 -0.9952 2.3708 19.5048
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 37.36119 0.22774 164.050 <2e-16 ***
## Depth 0.03351 0.02636 1.271 0.204
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.838 on 3157 degrees of freedom
## Multiple R-squared: 0.0005117, Adjusted R-squared: 0.0001951
## F-statistic: 1.616 on 1 and 3157 DF, p-value: 0.2037summary(lm(Temp.pav ~ factor(Method), data = DL))##
## Call:
## lm(formula = Temp.pav ~ factor(Method), data = DL)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.0672 -3.2737 -0.5672 2.6828 16.4328
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 37.0468 0.1328 278.883 <2e-16 ***
## factor(Method)NT 3.5204 0.1879 18.739 <2e-16 ***
## factor(Method)BB -1.7730 0.1879 -9.438 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.311 on 3156 degrees of freedom
## Multiple R-squared: 0.2068, Adjusted R-squared: 0.2063
## F-statistic: 411.4 on 2 and 3156 DF, p-value: < 2.2e-16summary(lm(Temp.pav ~ Temp.env, data = DL))##
## Call:
## lm(formula = Temp.pav ~ Temp.env, data = DL)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.1116 -2.6014 -0.1022 2.0891 15.1144
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 32.067132 0.137643 232.97 <2e-16 ***
## Temp.env 0.119006 0.002578 46.16 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.739 on 3157 degrees of freedom
## Multiple R-squared: 0.403, Adjusted R-squared: 0.4028
## F-statistic: 2131 on 1 and 3157 DF, p-value: < 2.2e-16# Multiple variables regression model
m1 <- summary(lm(Temp.pav ~ Hour + Depth + factor(Method)+ Temp.env, data = DL))
m1##
## Call:
## lm(formula = Temp.pav ~ Hour + Depth + factor(Method) + Temp.env,
## data = DL)
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.4733 -2.1042 -0.0814 1.9489 12.4180
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 30.275010 0.197993 152.909 <2e-16 ***
## Hour 0.098435 0.007700 12.784 <2e-16 ***
## Depth 0.033505 0.016058 2.087 0.037 *
## factor(Method)NT 3.520418 0.128461 27.405 <2e-16 ***
## factor(Method)BB -1.773029 0.128461 -13.802 <2e-16 ***
## Temp.env 0.112818 0.002089 54.001 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.948 on 3153 degrees of freedom
## Multiple R-squared: 0.6295, Adjusted R-squared: 0.6289
## F-statistic: 1071 on 5 and 3153 DF, p-value: < 2.2e-16## Xet them tuong tac giua bien Temp.env va Hour
m2 <- summary(lm(Temp.pav ~ Hour + Depth + factor(Method)+ Temp.env + Hour:Temp.env, data = DL))
m2##
## Call:
## lm(formula = Temp.pav ~ Hour + Depth + factor(Method) + Temp.env +
## Hour:Temp.env, data = DL)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.7633 -1.8774 -0.0705 1.7678 12.4405
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 33.1318927 0.2443428 135.596 <2e-16 ***
## Hour -0.1865645 0.0171791 -10.860 <2e-16 ***
## Depth 0.0335055 0.0152664 2.195 0.0283 *
## factor(Method)NT 3.5204179 0.1221314 28.825 <2e-16 ***
## factor(Method)BB -1.7730294 0.1221314 -14.517 <2e-16 ***
## Temp.env 0.0054494 0.0061826 0.881 0.3782
## Hour:Temp.env 0.0091412 0.0004985 18.338 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.802 on 3152 degrees of freedom
## Multiple R-squared: 0.6652, Adjusted R-squared: 0.6645
## F-statistic: 1044 on 6 and 3152 DF, p-value: < 2.2e-16m3 <- summary(lm(Temp.pav ~ Depth + factor(Method)+ Temp.env, data = DL))
m3##
## Call:
## lm(formula = Temp.pav ~ Depth + factor(Method) + Temp.env, data = DL)
##
## Residuals:
## Min 1Q Median 3Q Max
## -7.4132 -2.1510 -0.1583 2.0947 12.3105
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 31.216625 0.188455 165.645 <2e-16 ***
## Depth 0.033505 0.016466 2.035 0.042 *
## factor(Method)NT 3.520418 0.131727 26.725 <2e-16 ***
## factor(Method)BB -1.773029 0.131727 -13.460 <2e-16 ***
## Temp.env 0.119006 0.002084 57.105 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.023 on 3154 degrees of freedom
## Multiple R-squared: 0.6103, Adjusted R-squared: 0.6098
## F-statistic: 1235 on 4 and 3154 DF, p-value: < 2.2e-16# Find Parsimonious model - Tìm mô hình tối ưu
library(BMA)## Loading required package: leaps## Loading required package: robustbase##
## Attaching package: 'robustbase'## The following object is masked from 'package:survival':
##
## heart## Loading required package: inline## Loading required package: rrcov## Scalable Robust Estimators with High Breakdown Point (version 1.4-9)xvars <- cbind (DL$Hour, DL$Depth, DL$Method, DL$Temp.env)
y <- DL$Temp.pav
bma.search <- bicreg (xvars, y, strict = F, OR = 20)
summary(bma.search)##
## Call:
## bicreg(x = xvars, y = y, strict = F, OR = 20)
##
##
## 2 models were selected
## Best 2 models (cumulative posterior probability = 1 ):
##
## p!=0 EV SD model 1 model 2
## Intercept 100.0 32.879632 0.238622 3.290e+01 3.263e+01
## X1 100.0 0.098435 0.009425 9.843e-02 9.843e-02
## X2 7.1 0.002364 0.010043 . 3.351e-02
## X3 100.0 -0.886515 0.078615 -8.865e-01 -8.865e-01
## X4 100.0 0.112818 0.002557 1.128e-01 1.128e-01
##
## nVar 3 4
## r2 0.445 0.445
## BIC -1.833e+03 -1.828e+03
## post prob 0.929 0.071m4 <- summary(lm(Temp.pav ~ Hour + factor(Method)+ Temp.env, data = DL))
m4##
## Call:
## lm(formula = Temp.pav ~ Hour + factor(Method) + Temp.env, data = DL)
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.473 -2.114 -0.109 1.959 12.284
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 30.543054 0.150743 202.62 <2e-16 ***
## Hour 0.098435 0.007704 12.78 <2e-16 ***
## factor(Method)NT 3.520418 0.128529 27.39 <2e-16 ***
## factor(Method)BB -1.773029 0.128529 -13.79 <2e-16 ***
## Temp.env 0.112818 0.002090 53.97 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.949 on 3154 degrees of freedom
## Multiple R-squared: 0.629, Adjusted R-squared: 0.6285
## F-statistic: 1337 on 4 and 3154 DF, p-value: < 2.2e-16m5 <- summary(lm(Temp.pav ~ Hour + factor(Method)+ Temp.env + Hour:Temp.env, data = DL))
m5##
## Call:
## lm(formula = Temp.pav ~ Hour + factor(Method) + Temp.env + Hour:Temp.env,
## data = DL)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.7633 -1.8787 -0.0879 1.8055 12.3065
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 33.3999364 0.2117582 157.727 <2e-16 ***
## Hour -0.1865645 0.0171895 -10.853 <2e-16 ***
## factor(Method)NT 3.5204179 0.1222053 28.807 <2e-16 ***
## factor(Method)BB -1.7730294 0.1222053 -14.509 <2e-16 ***
## Temp.env 0.0054494 0.0061864 0.881 0.378
## Hour:Temp.env 0.0091412 0.0004988 18.327 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.804 on 3153 degrees of freedom
## Multiple R-squared: 0.6647, Adjusted R-squared: 0.6641
## F-statistic: 1250 on 5 and 3153 DF, p-value: < 2.2e-16