-
Notifications
You must be signed in to change notification settings - Fork 49
/
Module-5-Example-3.R
79 lines (51 loc) · 1.96 KB
/
Module-5-Example-3.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
# setwd("SET THE Working Director to THE PATH TO THIS DIRECTORY")
# emmeans:::convert_scripts()
# Load the smoker data set.
data<- read.csv("Datasets/smoking_SBP.csv")
# Print out the data and check what it is inside
print(data)
# View(data)
attach(data)
# Create dummy variables
data$g0 <- ifelse(data$group=='Current heavy smoker', 1, 0)
data$g1 <- ifelse(data$group=='Current light smoker', 1, 0)
data$g2 <- ifelse(data$group=='Former smoker', 1, 0)
data$g3 <- ifelse(data$group=='Never smoker', 1, 0)
print("Our Data after addomg Dummy Variables")
print(data)
# View(data)
# One-way ANOVA using lm() function
m2 <- lm(data$SBP~data$g0+data$g1+data$g2, data=data)
summary(m2)
m3 <- lm(data$SBP~data$g1+data$g2+data$g3, data=data)
summary(m3)
m4 <- lm(data$SBP~data$g0+data$g2+data$g3, data=data)
summary(m4)
# First run a normal one way anova
anova.model <- aov(data$SBP ~ group, data=data)
summary(anova.model)
# Install one time the "car" package.
# install.packages('car')
# Re-run ANOVA adjusting for Age
library(car)
# ANCOVA
# Now we run ANVOA with adjusting for age
Anova(lm(data$SBP ~ data$age + data$group), type=3)
# Type should be type 3 sum of squares.
# It defines the different types of sums of squares. Read more about it here
# https://www.r-bloggers.com/anova-%E2%80%93-type-iiiiii-ss-explained/
# Least square means
# install.packages('lsmeans') # one time installation only
library(lsmeans)
# Set the options for lsmeans
options(contrasts=c("contr.treatment", "contr.poly"))
# Calculate the lsmeans
lsmeans( lm(data$SBP~ data$group + data$age), pairwise~data$group, adjust="none")
# https://cran.r-project.org/web/packages/emmeans/vignettes/transition-from-lsmeans.html
# install.packages('emmeans')
library(emmeans)
my.model<-lm(SBP~group+age, data = data)
emm_options(contrasts=c("contr.treatment", "contr.poly"))
emmeans(my.model, specs = "group")
# or with pairwise comnparisions
emmeans(my.model, specs = "group", contr = "pairwise")