knitr::opts_chunk$set(echo = TRUE)
if(!require("tidyverse")) {install.packages("tidyverse"); library("tidyverse")}
if(!require("Hmisc")) {install.packages("Hmisc"); library("Hmisc")}
if(!require("sjlabelled")) {install.packages("sjlabelled"); library("sjlabelled")}
if(!require("DescTools")) {install.packages("DescTools"); library("DescTools")}
if(!require("GoodmanKruskal")) {install.packages("GoodmanKruskal"); library("GoodmanKruskal")}
if(!require("rcompanion")) {install.packages("rcompanion"); library("rcompanion")}
# clean up working environment
rm(list = ls())
# Generate data
mydatc6_case <- data.frame(Respondent = c(1:220),
Season = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4),
Type = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3))
# Set label
mydatc6_case$Season <- factor(mydatc6_case$Season,
levels = c(1:4),
labels = c("Spring",
"Summer",
"Fall",
"Winter"))
mydatc6_case$Type <- factor(mydatc6_case$Type,
levels = c(1:3),
labels = c("Milk",
"Dark",
"Yoghurt"))
print(mydatc6_case)
## Respondent Season Type
## 1 1 Spring Milk
## 2 2 Spring Milk
## 3 3 Spring Milk
## 4 4 Spring Milk
## 5 5 Spring Milk
## 6 6 Spring Milk
## 7 7 Spring Milk
## 8 8 Spring Milk
## 9 9 Spring Milk
## 10 10 Spring Milk
## 11 11 Spring Milk
## 12 12 Spring Milk
## 13 13 Spring Milk
## 14 14 Spring Milk
## 15 15 Spring Milk
## 16 16 Spring Milk
## 17 17 Spring Milk
## 18 18 Spring Milk
## 19 19 Spring Milk
## 20 20 Spring Dark
## 21 21 Spring Dark
## 22 22 Spring Dark
## 23 23 Spring Dark
## 24 24 Spring Dark
## 25 25 Spring Dark
## 26 26 Spring Dark
## 27 27 Spring Dark
## 28 28 Spring Yoghurt
## 29 29 Spring Yoghurt
## 30 30 Spring Yoghurt
## 31 31 Spring Yoghurt
## 32 32 Spring Yoghurt
## 33 33 Spring Yoghurt
## 34 34 Spring Yoghurt
## 35 35 Spring Yoghurt
## 36 36 Spring Yoghurt
## 37 37 Spring Yoghurt
## 38 38 Spring Yoghurt
## 39 39 Spring Yoghurt
## 40 40 Spring Yoghurt
## 41 41 Spring Yoghurt
## 42 42 Summer Milk
## 43 43 Summer Milk
## 44 44 Summer Milk
## 45 45 Summer Milk
## 46 46 Summer Milk
## 47 47 Summer Milk
## 48 48 Summer Milk
## 49 49 Summer Milk
## 50 50 Summer Milk
## 51 51 Summer Milk
## 52 52 Summer Milk
## 53 53 Summer Milk
## 54 54 Summer Milk
## 55 55 Summer Milk
## 56 56 Summer Dark
## 57 57 Summer Dark
## 58 58 Summer Dark
## 59 59 Summer Dark
## 60 60 Summer Dark
## 61 61 Summer Dark
## 62 62 Summer Yoghurt
## 63 63 Summer Yoghurt
## 64 64 Summer Yoghurt
## 65 65 Summer Yoghurt
## 66 66 Summer Yoghurt
## 67 67 Summer Yoghurt
## 68 68 Summer Yoghurt
## 69 69 Summer Yoghurt
## 70 70 Summer Yoghurt
## 71 71 Summer Yoghurt
## 72 72 Summer Yoghurt
## 73 73 Summer Yoghurt
## 74 74 Summer Yoghurt
## 75 75 Summer Yoghurt
## 76 76 Summer Yoghurt
## 77 77 Summer Yoghurt
## 78 78 Fall Milk
## 79 79 Fall Milk
## 80 80 Fall Milk
## 81 81 Fall Milk
## 82 82 Fall Milk
## 83 83 Fall Milk
## 84 84 Fall Milk
## 85 85 Fall Milk
## 86 86 Fall Milk
## 87 87 Fall Milk
## 88 88 Fall Milk
## 89 89 Fall Milk
## 90 90 Fall Milk
## 91 91 Fall Milk
## 92 92 Fall Milk
## 93 93 Fall Milk
## 94 94 Fall Milk
## 95 95 Fall Milk
## 96 96 Fall Milk
## 97 97 Fall Milk
## 98 98 Fall Milk
## 99 99 Fall Milk
## 100 100 Fall Milk
## 101 101 Fall Milk
## 102 102 Fall Milk
## 103 103 Fall Milk
## 104 104 Fall Milk
## 105 105 Fall Dark
## 106 106 Fall Dark
## 107 107 Fall Dark
## 108 108 Fall Dark
## 109 109 Fall Dark
## 110 110 Fall Dark
## 111 111 Fall Dark
## 112 112 Fall Dark
## 113 113 Fall Dark
## 114 114 Fall Dark
## 115 115 Fall Dark
## 116 116 Fall Dark
## 117 117 Fall Dark
## 118 118 Fall Dark
## 119 119 Fall Dark
## 120 120 Fall Dark
## 121 121 Fall Dark
## 122 122 Fall Dark
## 123 123 Fall Dark
## 124 124 Fall Dark
## 125 125 Fall Dark
## 126 126 Fall Dark
## 127 127 Fall Dark
## 128 128 Fall Dark
## 129 129 Fall Dark
## 130 130 Fall Dark
## 131 131 Fall Dark
## 132 132 Fall Yoghurt
## 133 133 Fall Yoghurt
## 134 134 Fall Yoghurt
## 135 135 Fall Yoghurt
## 136 136 Fall Yoghurt
## 137 137 Fall Yoghurt
## 138 138 Fall Yoghurt
## 139 139 Fall Yoghurt
## 140 140 Fall Yoghurt
## 141 141 Fall Yoghurt
## 142 142 Winter Milk
## 143 143 Winter Milk
## 144 144 Winter Milk
## 145 145 Winter Milk
## 146 146 Winter Milk
## 147 147 Winter Milk
## 148 148 Winter Milk
## 149 149 Winter Milk
## 150 150 Winter Milk
## 151 151 Winter Milk
## 152 152 Winter Milk
## 153 153 Winter Milk
## 154 154 Winter Milk
## 155 155 Winter Milk
## 156 156 Winter Milk
## 157 157 Winter Milk
## 158 158 Winter Milk
## 159 159 Winter Milk
## 160 160 Winter Milk
## 161 161 Winter Milk
## 162 162 Winter Milk
## 163 163 Winter Milk
## 164 164 Winter Milk
## 165 165 Winter Milk
## 166 166 Winter Milk
## 167 167 Winter Milk
## 168 168 Winter Milk
## 169 169 Winter Milk
## 170 170 Winter Milk
## 171 171 Winter Milk
## 172 172 Winter Milk
## 173 173 Winter Milk
## 174 174 Winter Milk
## 175 175 Winter Milk
## 176 176 Winter Milk
## 177 177 Winter Dark
## 178 178 Winter Dark
## 179 179 Winter Dark
## 180 180 Winter Dark
## 181 181 Winter Dark
## 182 182 Winter Dark
## 183 183 Winter Dark
## 184 184 Winter Dark
## 185 185 Winter Dark
## 186 186 Winter Dark
## 187 187 Winter Dark
## 188 188 Winter Dark
## 189 189 Winter Dark
## 190 190 Winter Dark
## 191 191 Winter Dark
## 192 192 Winter Dark
## 193 193 Winter Dark
## 194 194 Winter Dark
## 195 195 Winter Dark
## 196 196 Winter Dark
## 197 197 Winter Dark
## 198 198 Winter Dark
## 199 199 Winter Dark
## 200 200 Winter Dark
## 201 201 Winter Dark
## 202 202 Winter Dark
## 203 203 Winter Dark
## 204 204 Winter Dark
## 205 205 Winter Dark
## 206 206 Winter Dark
## 207 207 Winter Dark
## 208 208 Winter Dark
## 209 209 Winter Dark
## 210 210 Winter Dark
## 211 211 Winter Dark
## 212 212 Winter Dark
## 213 213 Winter Dark
## 214 214 Winter Dark
## 215 215 Winter Dark
## 216 216 Winter Yoghurt
## 217 217 Winter Yoghurt
## 218 218 Winter Yoghurt
## 219 219 Winter Yoghurt
## 220 220 Winter Yoghurt
Crosstable_Case_Chpt.6 <- table(mydatc6_case$Season, mydatc6_case$Type)
print(Crosstable_Case_Chpt.6);
##
## Milk Dark Yoghurt
## Spring 19 8 14
## Summer 14 6 16
## Fall 27 27 10
## Winter 35 39 5
margin.table(Crosstable_Case_Chpt.6, 1);
##
## Spring Summer Fall Winter
## 41 36 64 79
margin.table(Crosstable_Case_Chpt.6, 2)
##
## Milk Dark Yoghurt
## 95 80 45
addmargins(Crosstable_Case_Chpt.6, margin = 1);
##
## Milk Dark Yoghurt
## Spring 19 8 14
## Summer 14 6 16
## Fall 27 27 10
## Winter 35 39 5
## Sum 95 80 45
addmargins(Crosstable_Case_Chpt.6, margin = 2)
##
## Milk Dark Yoghurt Sum
## Spring 19 8 14 41
## Summer 14 6 16 36
## Fall 27 27 10 64
## Winter 35 39 5 79
# Figure 6.10 - Contingency analysis in SPSS (=R): Cross table
prop.table(Crosstable_Case_Chpt.6)
##
## Milk Dark Yoghurt
## Spring 0.08636364 0.03636364 0.06363636
## Summer 0.06363636 0.02727273 0.07272727
## Fall 0.12272727 0.12272727 0.04545455
## Winter 0.15909091 0.17727273 0.02272727
prop.table(Crosstable_Case_Chpt.6, m = 1)
##
## Milk Dark Yoghurt
## Spring 0.46341463 0.19512195 0.34146341
## Summer 0.38888889 0.16666667 0.44444444
## Fall 0.42187500 0.42187500 0.15625000
## Winter 0.44303797 0.49367089 0.06329114
prop.table(Crosstable_Case_Chpt.6, m = 2)
##
## Milk Dark Yoghurt
## Spring 0.2000000 0.1000000 0.3111111
## Summer 0.1473684 0.0750000 0.3555556
## Fall 0.2842105 0.3375000 0.2222222
## Winter 0.3684211 0.4875000 0.1111111
summary(Crosstable_Case_Chpt.6);
## Number of cases in table: 220
## Number of factors: 2
## Test for independence of all factors:
## Chisq = 33.92, df = 6, p-value = 6.965e-06
chi <- chisq.test(Crosstable_Case_Chpt.6)
print(chi);
##
## Pearson's Chi-squared test
##
## data: Crosstable_Case_Chpt.6
## X-squared = 33.922, df = 6, p-value = 6.965e-06
chi$expected;
##
## Milk Dark Yoghurt
## Spring 17.70455 14.90909 8.386364
## Summer 15.54545 13.09091 7.363636
## Fall 27.63636 23.27273 13.090909
## Winter 34.11364 28.72727 16.159091
chi$residuals
##
## Milk Dark Yoghurt
## Spring 0.3078788 -1.7893501 1.9384626
## Summer -0.3919715 -1.9598237 3.1826197
## Fall -0.1210500 0.7726228 -0.8542821
## Winter 0.1517567 1.9166312 -2.7760057
chi$p.value
## [1] 6.964982e-06
phi <- sqrt(as.numeric(unlist(chi$statistic))/sum(Crosstable_Case_Chpt.6))
Contingency_Coefficient <- sqrt(as.numeric(unlist(chi$statistic))/(sum(Crosstable_Case_Chpt.6) + as.numeric(unlist(chi$statistic))))
print(phi);
## [1] 0.3926711
rcompanion::cramerV(Crosstable_Case_Chpt.6);
## Cramer V
## 0.2777
print(Contingency_Coefficient)
## [1] 0.3655023
barplot(Crosstable_Case_Chpt.6,
beside = TRUE,
main = "Barplot of the cross tabulation",
axes = TRUE,
xlab = "Type",
ylab = "Volume",
col = c("red","blue", "orange", "grey"),
legend.text = c("Spring","Summer", "Fall", "Winter"),
names.arg = c("Milk","Dark", "Yoghurt"),
cex.main = 1,
cex.lab = 1,
cex.axis = 1)
# Lambda (symmetric)
paste("Lambda (symmetric)");
## [1] "Lambda (symmetric)"
DescTools::Lambda(mydatc6_case$Season, mydatc6_case$Type, direction = "symmetric", conf.level = 0.95);
## lambda lwr.ci upr.ci
## 0.06390977 0.00000000 0.15355226
# Lambda (Season dependent)
paste("Lambda (Season dependent)");
## [1] "Lambda (Season dependent)"
DescTools::Lambda(mydatc6_case$Season, mydatc6_case$Type, direction = "row", conf.level = 0.95);
## lambda lwr.ci upr.ci
## 0.07801418 0.01684950 0.13917887
# Lambda (Type dependent)
paste("Lambda (Type dependent)");
## [1] "Lambda (Type dependent)"
DescTools::Lambda(mydatc6_case$Season, mydatc6_case$Type, direction = "column", conf.level = 0.95)
## lambda lwr.ci upr.ci
## 0.0480000 0.0000000 0.2040175
# GoodmanKruskal (tauxy = Type dependent; tauyx = Season dependent)
paste("GoodmanKruskal (tauxy = Type dependent; tauyx = Season dependent)");
## [1] "GoodmanKruskal (tauxy = Type dependent; tauyx = Season dependent)"
GoodmanKruskal::GKtau(mydatc6_case$Season, mydatc6_case$Type)
## xName yName Nx Ny tauxy tauyx
## 1 mydatc6_case$Season mydatc6_case$Type 4 3 0.063 0.051