# Loading the library
library(plyr)
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:plyr':
##
## arrange, count, desc, failwith, id, mutate, rename, summarise,
## summarize
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(tidyverse)
## -- Attaching packages --------------------------- tidyverse 1.3.0 --
## v ggplot2 3.3.2 v purrr 0.3.4
## v tibble 3.0.3 v stringr 1.4.0
## v tidyr 1.1.2 v forcats 0.5.0
## v readr 1.3.1
## -- Conflicts ------------------------------ tidyverse_conflicts() --
## x dplyr::arrange() masks plyr::arrange()
## x purrr::compact() masks plyr::compact()
## x dplyr::count() masks plyr::count()
## x dplyr::failwith() masks plyr::failwith()
## x dplyr::filter() masks stats::filter()
## x dplyr::id() masks plyr::id()
## x dplyr::lag() masks stats::lag()
## x dplyr::mutate() masks plyr::mutate()
## x dplyr::rename() masks plyr::rename()
## x dplyr::summarise() masks plyr::summarise()
## x dplyr::summarize() masks plyr::summarize()
library(tidyr)
library(sjPlot)
## Warning: package 'sjPlot' was built under R version 4.0.3
library(sjmisc)
## Warning: package 'sjmisc' was built under R version 4.0.3
##
## Attaching package: 'sjmisc'
## The following object is masked from 'package:purrr':
##
## is_empty
## The following object is masked from 'package:tidyr':
##
## replace_na
## The following object is masked from 'package:tibble':
##
## add_case
library(sjlabelled)
## Warning: package 'sjlabelled' was built under R version 4.0.3
##
## Attaching package: 'sjlabelled'
## The following objects are masked from 'package:sjmisc':
##
## to_character, to_factor, to_label, to_numeric
## The following object is masked from 'package:forcats':
##
## as_factor
## The following object is masked from 'package:dplyr':
##
## as_label
library(outreg)
## Warning: package 'outreg' was built under R version 4.0.3
# Loading the data
library(readxl)
data <- read_excel("C:/Users/ramya.emandi/Desktop/Econ Policy/PROJECTS/MSME Study/for R_labelled.xlsx")
colnames(data)[7] <- "Age"
data <- data %>% mutate(`MW - Male` = as.numeric(`MW - Male`))
## Warning: Problem with `mutate()` input `MW - Male`.
## i NAs introduced by coercion
## i Input `MW - Male` is `as.numeric(`MW - Male`)`.
## Warning in mask$eval_all_mutate(dots[[i]]): NAs introduced by coercion
data <- data %>% mutate(`MW - Female` = as.numeric(`MW - Female`))
## Warning: Problem with `mutate()` input `MW - Female`.
## i NAs introduced by coercion
## i Input `MW - Female` is `as.numeric(`MW - Female`)`.
## Warning: NAs introduced by coercion
data <- data %>% mutate(prc_loss = as.numeric(prc_loss))
## Warning: Problem with `mutate()` input `prc_loss`.
## i NAs introduced by coercion
## i Input `prc_loss` is `as.numeric(prc_loss)`.
## Warning: NAs introduced by coercion
data <- data %>% mutate(Sector_Company = as.factor(Sector_Company))
levels(data$Sector_Company) <- list(AS = "Agriculture Sector", MS = "Manufacturing Sector", SS = "Service Sector")
data <- data %>% mutate(Type_Company = as.factor(Type_Company))
levels(data$Type_Company) <- list(CLG = "Company Limited by Guarantee",CLS_Pvt = "Company Limited by Shares: Private Company", CLS_Pub = "Company Limited by Shares: Public Company", UL = "Unlimited Company")
data <- data %>% mutate(Size_Company = as.factor(Size_Company))
data$Size_Company <- revalue(data$Size_Company, c("Large (Investment in Plant and Machinery or Equipment: more than Rs.50 crore and Annual Turnover ; more than Rs. 250 crore)" = "large",
"Medium (Investment in Plant and Machinery or Equipment: Not more than Rs.50 crore and Annual Turnover ; not more than Rs. 250 crore)" = "Medium",
"Micro (Investment in Plant and Machinery or Equipment: Not more than Rs.1 crore and Annual Turnover ; not more than Rs. 5 crore)" = "Micro",
"Others" ="Others",
"Self Employed" = "Self",
"Small (Investment in Plant and Machinery or Equipment: Not more than Rs.10 crore and Annual Turnover ; not more than Rs. 50 crore)" = "Small"))
data$Size_Company <- factor(data$Size_Company, ordered = TRUE, levels = c("Self", "Micro", "Small", "Medium", "Large"))
data <- data %>% mutate(impact_severity = as.factor(impact_severity))
data$impact_severity <- revalue(data$impact_severity, c("Marginal, Will be back to normal operations and business in next few months" = "Marginal",
"Moderate, WIll be back to normal operations and business in next few years" = "Moderate",
"None, Back to normal operations and business" = "None",
"Others" = "Others",
"Severe, Cannot recover/Shutting down" = "Severe"))
#imp_ord <- ordered(c("None", "Marginal", "Moderate", "Severe", "Others"))
data$Impacts <- factor(data$impact_severity, ordered = TRUE, levels = c("None", "Marginal", "Moderate", "Severe", "Others"))
data$impact_severity <- factor(data$impact_severity, ordered = TRUE, levels = c("None", "Marginal", "Moderate", "Severe"))
levels(data$impact_severity)
## [1] "None" "Marginal" "Moderate" "Severe"
data <- data %>% mutate(Caste = as.factor(Caste))
data <- data %>% mutate(duration_skilltrain = as.factor(duration_skilltrain))
data <- data %>% mutate(R_U = as.factor(R_U))
data <- data %>% mutate(Gender = as.factor(Gender))
data <- data %>% mutate(Diff_abled = as.factor(Diff_abled))
data <- data %>% mutate(industrial_associatn = as.factor(industrial_associatn))
data <- data %>% mutate(SEZ = as.factor(SEZ))
data <- data %>% mutate(aware_govt_schemes = as.factor(aware_govt_schemes))
data <- data %>% mutate(DYKaccess_govt_schemes = as.factor(DYKaccess_govt_schemes))
data <- data %>% mutate(accessing_govt_schemes = as.factor(accessing_govt_schemes))
data <- data %>% mutate(interest_webinars = as.factor(interest_webinars))
data <- data %>% mutate(interest_skilltrain = as.factor(interest_skilltrain))
data <- data %>% mutate(MW_left_covid = as.factor(MW_left_covid))
View(data)
#Multinomial Logistic Regression
#setting the baseline
#data$Sector_Company <- relevel(data$Sector_Company)
#ml <- multinom(formula = Sector_Company ~ ILBC + ILLC, data = data)
#summary(ml)
#Descriptve GRAPHS
###########################################################
continuous <- c("Age", "Total_Employees", "Migrant_Employees", "Vacancies_Employees", "prc_loss")
#summary(data$impact_severity, data$Sector_Company, data$Type_Company, data$Size_Company, data$Age, data$R_U)
summary(data[continuous])
## Age Total_Employees Migrant_Employees Vacancies_Employees
## Min. : 1.00 Min. : 0.00 Min. : 0.000 Min. : 0.000
## 1st Qu.: 6.00 1st Qu.: 5.00 1st Qu.: 0.000 1st Qu.: 0.000
## Median :12.00 Median : 9.00 Median : 0.000 Median : 0.000
## Mean :14.62 Mean : 21.41 Mean : 3.585 Mean : 2.462
## 3rd Qu.:20.00 3rd Qu.: 18.00 3rd Qu.: 1.000 3rd Qu.: 0.000
## Max. :86.00 Max. :1500.00 Max. :500.000 Max. :460.000
##
## prc_loss
## Min. : 0.00
## 1st Qu.: 40.00
## Median : 50.00
## Mean : 50.31
## 3rd Qu.: 60.00
## Max. :560.00
## NA's :1
library(pastecs)
## Warning: package 'pastecs' was built under R version 4.0.3
##
## Attaching package: 'pastecs'
## The following object is masked from 'package:tidyr':
##
## extract
## The following objects are masked from 'package:dplyr':
##
## first, last
res <- (stat.desc(data[continuous]))
res <- round(res, 2)
View(res)
#Discrete <- c("MW_left_covid", "R_U", "Gender", "Caste", "Diff_abled", "industrial_associatn", "Sector_Company", "Size_Company", "Type_Company", "SEZ")
#summary(data[discrete])
plot <- ggplot(data=data, aes(x=Age, y=Total_Employees, col=impact_severity))
plot <- plot + geom_point(aes(size = 1.5))
plot <- plot + xlab("Total Employees") + ylab("% loss") +
scale_color_discrete(name = "COVID Impact") + ylim(c(0,200))
plot
## Warning: Removed 12 rows containing missing values (geom_point).
ratiom_t <- sum(data$Migrant_Employees)/ sum(data$Total_Employees)
ratiom_t
## [1] 0.1674665
mean(data$prc_loss)
## [1] NA
library("ggpubr")
## Warning: package 'ggpubr' was built under R version 4.0.3
## 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: 'ggpubr'
## The following object is masked from 'package:plyr':
##
## mutate
ggboxplot(data, x = "impact_severity", y = "prc_loss",
color = "impact_severity", ylim = c(0,100),
ylab = "% loss", xlab = "Impact")
## Warning: Removed 1 rows containing non-finite values (stat_boxplot).
ggboxplot(data, x = "impact_severity", y = "Age",
color = "impact_severity", ylim = c(0,30),
ylab = "Age of the Company", xlab = "Impact")
ggboxplot(data, x = "impact_severity", y = "Migrant_Employees",
color = "impact_severity", ylim = c(0,5),
ylab = "Employees left", xlab = "Impact")
##########################################
ggplot(data, aes(y = Age)) +
coord_cartesian(ylim = c(0, 80)) +
geom_boxplot(size = .75) +
facet_grid(R_U ~ Sector_Company, margins = TRUE) +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1))
ggplot(data, aes(y = Migrant_Employees)) +
coord_cartesian(ylim = c(0, 20)) +
geom_boxplot(size = .75) +
facet_grid(R_U ~ Sector_Company, margins = TRUE) +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1))
ggplot(data, aes(y = Total_Employees)) +
coord_cartesian(ylim = c(0, 30)) +
geom_boxplot(size = .75) +
facet_grid(R_U ~ Sector_Company, margins = TRUE) +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1))
#migrants employed versus migrants that left
data <- data %>% mutate(return_MW = as.factor(return_MW))
ggplot(data, aes(x = MW_left_covid, y = Migrant_Employees)) +
coord_cartesian(ylim = c(0, 60)) +
geom_boxplot(size = .75) +
facet_grid(R_U ~ Sector_Company, margins = TRUE) +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1))
#linear regression
#########################################################
lin <- lm(formula = prc_loss ~ MW_left_covid + Age + Sector_Company + Type_Company, data = data)
summary(lin)
##
## Call:
## lm(formula = prc_loss ~ MW_left_covid + Age + Sector_Company +
## Type_Company, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -55.56 -14.09 -0.62 12.57 521.07
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 43.91274 7.94613 5.526 4.24e-08 ***
## MW_left_covidYes 5.13104 1.94204 2.642 0.00838 **
## Age -0.08538 0.07707 -1.108 0.26822
## Sector_CompanyMS 15.45095 2.66075 5.807 8.71e-09 ***
## Sector_CompanySS 15.23879 3.45181 4.415 1.13e-05 ***
## Type_CompanyCLS_Pvt -4.89880 7.53400 -0.650 0.51571
## Type_CompanyCLS_Pub -35.46626 13.48436 -2.630 0.00867 **
## Type_CompanyUL -8.20826 7.50288 -1.094 0.27423
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 27.61 on 936 degrees of freedom
## (10 observations deleted due to missingness)
## Multiple R-squared: 0.05133, Adjusted R-squared: 0.04424
## F-statistic: 7.236 on 7 and 936 DF, p-value: 1.846e-08
#ANOVA
#hist(data$prc_loss, breaks = c(45,50,55,60,65,70,75,80,85,90,95,100))
two.way <- aov(prc_loss ~ MW_left_covid + Sector_Company + Type_Company, data = data)
#outreg(two.way)
t <- summary(two.way)
#tab_model(t)
plot(two.way)
tukey.two.way<-TukeyHSD(two.way)
#table(tukey.two.way)
ggplot(data, aes(x = MW_left_covid, y = prc_loss)) +
coord_cartesian(ylim = c(0, 30)) +
geom_boxplot(size = .75) +
facet_wrap(~ Sector_Company) +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1))
## Warning: Removed 1 rows containing non-finite values (stat_boxplot).
# Logistic Regression
##########################################################
library(aod)
## Warning: package 'aod' was built under R version 4.0.3
library(ggplot2)
#data <- rename(data, c())
#data <- data %>% mutate(MW_left_covid = as.numeric(MW_left_covid))
#data$MW_left_covid <- data$MW_left_covid - 1
lr <- glm(formula = MW_left_covid ~ Total_Employees + Impacts + Age + R_U + Sector_Company, data = data, family = "binomial")
summary(lr)
##
## Call:
## glm(formula = MW_left_covid ~ Total_Employees + Impacts + Age +
## R_U + Sector_Company, family = "binomial", data = data)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.1883 -0.9021 -0.7412 1.2858 2.4035
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.001744 0.407842 -0.004 0.996588
## Total_Employees 0.007061 0.002205 3.202 0.001364 **
## Impacts.L 0.310312 0.583268 0.532 0.594710
## Impacts.Q -1.189392 0.498480 -2.386 0.017031 *
## Impacts.C -1.167942 0.411240 -2.840 0.004511 **
## Impacts^4 -0.634732 0.263479 -2.409 0.015995 *
## Age 0.021840 0.006072 3.597 0.000322 ***
## R_UUrban -0.790529 0.308881 -2.559 0.010487 *
## Sector_CompanyMS -0.308957 0.196421 -1.573 0.115735
## Sector_CompanySS -1.461629 0.324713 -4.501 6.75e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1226.6 on 953 degrees of freedom
## Residual deviance: 1139.7 on 944 degrees of freedom
## AIC: 1159.7
##
## Number of Fisher Scoring iterations: 5
tab_model(lr)
|
|
MW_left_covid
|
|
Predictors
|
Odds Ratios
|
CI
|
p
|
|
(Intercept)
|
1.00
|
0.44 – 2.21
|
0.997
|
|
Total_Employees
|
1.01
|
1.00 – 1.01
|
0.001
|
|
Impacts [linear]
|
1.36
|
0.37 – 3.99
|
0.595
|
|
Impacts [quadratic]
|
0.30
|
0.10 – 0.75
|
0.017
|
|
Impacts [cubic]
|
0.31
|
0.13 – 0.68
|
0.005
|
|
Impacts [4th degree]
|
0.53
|
0.31 – 0.88
|
0.016
|
|
Age
|
1.02
|
1.01 – 1.03
|
<0.001
|
|
R_U [Urban]
|
0.45
|
0.25 – 0.83
|
0.010
|
|
Sector_Company [MS]
|
0.73
|
0.50 – 1.08
|
0.116
|
|
Sector_Company [SS]
|
0.23
|
0.12 – 0.43
|
<0.001
|
|
Observations
|
954
|
|
R2 Tjur
|
0.088
|
o <- outreg(lr, pv = TRUE, tv = TRUE, starred = 'se')
View(o)
# wald test for overall signifince of categorical variable
wald.test(b = coef(lr), Sigma = vcov(lr), Terms = 7:8)
## Wald test:
## ----------
##
## Chi-squared test:
## X2 = 19.7, df = 2, P(> X2) = 5.3e-05
# odds ratios only
exp(coef(lr))
## (Intercept) Total_Employees Impacts.L Impacts.Q
## 0.9982574 1.0070856 1.3638506 0.3044063
## Impacts.C Impacts^4 Age R_UUrban
## 0.3110065 0.5300773 1.0220807 0.4536046
## Sector_CompanyMS Sector_CompanySS
## 0.7342125 0.2318582
# Confidence interval
exp(cbind(OR = coef(lr), confint(lr)))
## Waiting for profiling to be done...
## OR 2.5 % 97.5 %
## (Intercept) 0.9982574 0.44430029 2.2142513
## Total_Employees 1.0070856 1.00309163 1.0118078
## Impacts.L 1.3638506 0.36891619 3.9888711
## Impacts.Q 0.3044063 0.09936803 0.7517693
## Impacts.C 0.3110065 0.13157061 0.6752054
## Impacts^4 0.5300773 0.30952933 0.8766867
## Age 1.0220807 1.01006844 1.0344404
## R_UUrban 0.4536046 0.24610772 0.8312851
## Sector_CompanyMS 0.7342125 0.49989257 1.0805842
## Sector_CompanySS 0.2318582 0.11990136 0.4301909
contrasts(data$MW_left_covid)
## Yes
## No 0
## Yes 1
contrasts(data$impact_severity)
## .L .Q .C
## [1,] -0.6708204 0.5 -0.2236068
## [2,] -0.2236068 -0.5 0.6708204
## [3,] 0.2236068 -0.5 -0.6708204
## [4,] 0.6708204 0.5 0.2236068
# Predicted probabilities of migrants leaving at each value of Impact severity, holding others at their means
lr1 <- glm(formula = MW_left_covid ~ Total_Employees + Impacts + Age, data = data, family = "binomial")
newdata1 <- with(data, data.frame(Age = mean(Age), Total_Employees = mean(Total_Employees), Impacts = factor(c("None", "Marginal", "Moderate", "Severe", "Others"))))
newdata1$P <- predict(lr1, newdata = newdata1, type = "response")
newdata1
## Age Total_Employees Impacts P
## 1 14.6195 21.40671 None 0.2446020
## 2 14.6195 21.40671 Marginal 0.3160059
## 3 14.6195 21.40671 Moderate 0.3666538
## 4 14.6195 21.40671 Severe 0.6846303
## 5 14.6195 21.40671 Others 0.1819949
#Impact Severity
newdata2 <- with(data, data.frame(Age = rep(seq(from = 6, to = 90, length.out = 15), 5),
Total_Employees = mean(Total_Employees), Impacts = factor(rep(c("None", "Marginal", "Moderate", "Severe", "Others"), each = 15))))
newdata2
## Age Total_Employees Impacts
## 1 6 21.40671 None
## 2 12 21.40671 None
## 3 18 21.40671 None
## 4 24 21.40671 None
## 5 30 21.40671 None
## 6 36 21.40671 None
## 7 42 21.40671 None
## 8 48 21.40671 None
## 9 54 21.40671 None
## 10 60 21.40671 None
## 11 66 21.40671 None
## 12 72 21.40671 None
## 13 78 21.40671 None
## 14 84 21.40671 None
## 15 90 21.40671 None
## 16 6 21.40671 Marginal
## 17 12 21.40671 Marginal
## 18 18 21.40671 Marginal
## 19 24 21.40671 Marginal
## 20 30 21.40671 Marginal
## 21 36 21.40671 Marginal
## 22 42 21.40671 Marginal
## 23 48 21.40671 Marginal
## 24 54 21.40671 Marginal
## 25 60 21.40671 Marginal
## 26 66 21.40671 Marginal
## 27 72 21.40671 Marginal
## 28 78 21.40671 Marginal
## 29 84 21.40671 Marginal
## 30 90 21.40671 Marginal
## 31 6 21.40671 Moderate
## 32 12 21.40671 Moderate
## 33 18 21.40671 Moderate
## 34 24 21.40671 Moderate
## 35 30 21.40671 Moderate
## 36 36 21.40671 Moderate
## 37 42 21.40671 Moderate
## 38 48 21.40671 Moderate
## 39 54 21.40671 Moderate
## 40 60 21.40671 Moderate
## 41 66 21.40671 Moderate
## 42 72 21.40671 Moderate
## 43 78 21.40671 Moderate
## 44 84 21.40671 Moderate
## 45 90 21.40671 Moderate
## 46 6 21.40671 Severe
## 47 12 21.40671 Severe
## 48 18 21.40671 Severe
## 49 24 21.40671 Severe
## 50 30 21.40671 Severe
## 51 36 21.40671 Severe
## 52 42 21.40671 Severe
## 53 48 21.40671 Severe
## 54 54 21.40671 Severe
## 55 60 21.40671 Severe
## 56 66 21.40671 Severe
## 57 72 21.40671 Severe
## 58 78 21.40671 Severe
## 59 84 21.40671 Severe
## 60 90 21.40671 Severe
## 61 6 21.40671 Others
## 62 12 21.40671 Others
## 63 18 21.40671 Others
## 64 24 21.40671 Others
## 65 30 21.40671 Others
## 66 36 21.40671 Others
## 67 42 21.40671 Others
## 68 48 21.40671 Others
## 69 54 21.40671 Others
## 70 60 21.40671 Others
## 71 66 21.40671 Others
## 72 72 21.40671 Others
## 73 78 21.40671 Others
## 74 84 21.40671 Others
## 75 90 21.40671 Others
newdata3 <- cbind(newdata2, predict(lr1, newdata = newdata2, type = "link", se = TRUE))
newdata3 <- within(newdata3, {
PredictedProb <- plogis(fit)
LL <- plogis(fit - (1.96 * se.fit))
UL <- plogis(fit + (1.96 * se.fit))
})
newdata3
## Age Total_Employees Impacts fit se.fit residual.scale UL
## 1 6 21.40671 None -1.33916124 0.4064859 1 0.3676181
## 2 12 21.40671 None -1.19190292 0.4021020 1 0.4004038
## 3 18 21.40671 None -1.04464460 0.4008306 1 0.4356054
## 4 24 21.40671 None -0.89738628 0.4027011 1 0.4730033
## 5 30 21.40671 None -0.75012797 0.4076704 1 0.5122241
## 6 36 21.40671 None -0.60286965 0.4156274 1 0.5527431
## 7 42 21.40671 None -0.45561133 0.4264047 1 0.5939073
## 8 48 21.40671 None -0.30835301 0.4397951 1 0.6349809
## 9 54 21.40671 None -0.16109470 0.4555682 1 0.6752043
## 10 60 21.40671 None -0.01383638 0.4734860 1 0.7138581
## 11 66 21.40671 None 0.13342194 0.4933149 1 0.7503199
## 12 72 21.40671 None 0.28068026 0.5148340 1 0.7841057
## 13 78 21.40671 None 0.42793857 0.5378405 1 0.8148905
## 14 84 21.40671 None 0.57519689 0.5621519 1 0.8425088
## 15 90 21.40671 None 0.72245521 0.5876062 1 0.8669393
## 16 6 21.40671 Marginal -0.98373709 0.1154416 1 0.3191955
## 17 12 21.40671 Marginal -0.83647877 0.1030585 1 0.3464945
## 18 18 21.40671 Marginal -0.68922045 0.1021727 1 0.3801381
## 19 24 21.40671 Marginal -0.54196213 0.1130548 1 0.4205844
## 20 30 21.40671 Marginal -0.39470381 0.1328437 1 0.4664679
## 21 36 21.40671 Marginal -0.24744550 0.1582323 1 0.5156673
## 22 42 21.40671 Marginal -0.10018718 0.1869528 1 0.5661697
## 23 48 21.40671 Marginal 0.04707114 0.2176906 1 0.6162697
## 24 54 21.40671 Marginal 0.19432946 0.2497018 1 0.6645740
## 25 60 21.40671 Marginal 0.34158777 0.2825538 1 0.7100019
## 26 66 21.40671 Marginal 0.48884609 0.3159846 1 0.7517889
## 27 72 21.40671 Marginal 0.63610441 0.3498283 1 0.7894757
## 28 78 21.40671 Marginal 0.78336273 0.3839756 1 0.8228759
## 29 84 21.40671 Marginal 0.93062104 0.4183524 1 0.8520274
## 30 90 21.40671 Marginal 1.07787936 0.4529062 1 0.8771351
## 31 6 21.40671 Moderate -0.75814777 0.1185617 1 0.3715047
## 32 12 21.40671 Moderate -0.61088946 0.1069317 1 0.4009992
## 33 18 21.40671 Moderate -0.46363114 0.1064698 1 0.4366054
## 34 24 21.40671 Moderate -0.31637282 0.1173080 1 0.4784012
## 35 30 21.40671 Moderate -0.16911450 0.1367862 1 0.5247264
## 36 36 21.40671 Moderate -0.02185619 0.1618137 1 0.5732928
## 37 42 21.40671 Moderate 0.12540213 0.1902125 1 0.6220406
## 38 48 21.40671 Moderate 0.27266045 0.2206850 1 0.6693403
## 39 54 21.40671 Moderate 0.41991877 0.2524814 1 0.7139778
## 40 60 21.40671 Moderate 0.56717709 0.2851592 1 0.7551164
## 41 66 21.40671 Moderate 0.71443540 0.3184472 1 0.7922583
## 42 72 21.40671 Moderate 0.86169372 0.3521724 1 0.8251954
## 43 78 21.40671 Moderate 1.00895204 0.3862202 1 0.8539525
## 44 84 21.40671 Moderate 1.15621036 0.4205124 1 0.8787254
## 45 90 21.40671 Moderate 1.30346867 0.4549936 1 0.8998227
## 46 6 21.40671 Severe 0.56358448 0.4416549 1 0.8067810
## 47 12 21.40671 Severe 0.71084280 0.4396969 1 0.8281609
## 48 18 21.40671 Severe 0.85810112 0.4406036 1 0.8483457
## 49 24 21.40671 Severe 1.00535943 0.4443576 1 0.8671856
## 50 30 21.40671 Severe 1.15261775 0.4508877 1 0.8845619
## 51 36 21.40671 Severe 1.29987607 0.4600758 1 0.9003953
## 52 42 21.40671 Severe 1.44713439 0.4717665 1 0.9146512
## 53 48 21.40671 Severe 1.59439270 0.4857793 1 0.9273394
## 54 54 21.40671 Severe 1.74165102 0.5019196 1 0.9385097
## 55 60 21.40671 Severe 1.88890934 0.5199894 1 0.9482448
## 56 66 21.40671 Severe 2.03616766 0.5397949 1 0.9566514
## 57 72 21.40671 Severe 2.18342597 0.5611525 1 0.9638509
## 58 78 21.40671 Severe 2.33068429 0.5838917 1 0.9699713
## 59 84 21.40671 Severe 2.47794261 0.6078577 1 0.9751406
## 60 90 21.40671 Severe 2.62520093 0.6329110 1 0.9794817
## 61 6 21.40671 Others -1.71443882 0.7911968 1 0.4591679
## 62 12 21.40671 Others -1.56718051 0.7883128 1 0.4944784
## 63 18 21.40671 Others -1.41992219 0.7870233 1 0.5306225
## 64 24 21.40671 Others -1.27266387 0.7873362 1 0.5672194
## 65 30 21.40671 Others -1.12540555 0.7892495 1 0.6038478
## 66 36 21.40671 Others -0.97814724 0.7927517 1 0.6400650
## 67 42 21.40671 Others -0.83088892 0.7978218 1 0.6754286
## 68 48 21.40671 Others -0.68363060 0.8044302 1 0.7095197
## 69 54 21.40671 Others -0.53637228 0.8125393 1 0.7419646
## 70 60 21.40671 Others -0.38911396 0.8221048 1 0.7724525
## 71 66 21.40671 Others -0.24185565 0.8330764 1 0.8007477
## 72 72 21.40671 Others -0.09459733 0.8453995 1 0.8266954
## 73 78 21.40671 Others 0.05266099 0.8590158 1 0.8502206
## 74 84 21.40671 Others 0.19991931 0.8738650 1 0.8713216
## 75 90 21.40671 Others 0.34717762 0.8898853 1 0.8900599
## LL PredictedProb
## 1 0.10565868 0.2076480
## 2 0.12131641 0.2329188
## 3 0.13820581 0.2602548
## 4 0.15621289 0.2895879
## 5 0.17520733 0.3207934
## 6 0.19505392 0.3536874
## 7 0.21562194 0.3880275
## 8 0.23679148 0.4235168
## 9 0.25845612 0.4598132
## 10 0.28052297 0.4965410
## 11 0.30291075 0.5333061
## 12 0.32554717 0.5697130
## 13 0.34836612 0.6053813
## 14 0.37130537 0.6399615
## 15 0.39430475 0.6731474
## 16 0.22970059 0.2721509
## 17 0.26144467 0.3022769
## 18 0.29121737 0.3342065
## 19 0.31787611 0.3677313
## 20 0.34184624 0.4025855
## 21 0.36410742 0.4384523
## 22 0.38541780 0.4749741
## 23 0.40622278 0.5117656
## 24 0.42675919 0.5484301
## 25 0.44714379 0.5845762
## 26 0.46742526 0.6198346
## 27 0.48761278 0.6538723
## 28 0.50769201 0.6864044
## 29 0.52763442 0.7172013
## 30 0.54740302 0.7460925
## 31 0.27080766 0.3190485
## 32 0.30566273 0.3518563
## 33 0.33797936 0.3861248
## 34 0.36672407 0.4215600
## 35 0.39240469 0.4578218
## 36 0.41604974 0.4945362
## 37 0.43846001 0.5313095
## 38 0.46011440 0.5677459
## 39 0.48127258 0.6034638
## 40 0.50206626 0.6381115
## 41 0.52255442 0.6713805
## 42 0.54275431 0.7030144
## 43 0.56265897 0.7328150
## 44 0.58224740 0.7606434
## 45 0.60149091 0.7864182
## 46 0.42505418 0.6372815
## 47 0.46233073 0.6705874
## 48 0.49862950 0.7022638
## 49 0.53355412 0.7321110
## 50 0.56681738 0.7599887
## 51 0.59823769 0.7858141
## 52 0.62772562 0.8095570
## 53 0.65526537 0.8312332
## 54 0.68089517 0.8508967
## 55 0.70468955 0.8686311
## 56 0.72674488 0.8845425
## 57 0.74716844 0.8987513
## 58 0.76607087 0.9113866
## 59 0.78356124 0.9225810
## 60 0.79974405 0.9324660
## 61 0.03678507 0.1525889
## 62 0.04260381 0.1726187
## 63 0.04914960 0.1946738
## 64 0.05647404 0.2188016
## 65 0.06462570 0.2450100
## 66 0.07364914 0.2732596
## 67 0.08358390 0.3034571
## 68 0.09446353 0.3354515
## 69 0.10631478 0.3690319
## 70 0.11915681 0.4039306
## 71 0.13300051 0.4398291
## 72 0.14784799 0.4763683
## 73 0.16369204 0.5131622
## 74 0.18051591 0.5498140
## 75 0.19829299 0.5859330
ggplot(newdata3, aes(x = Age, y = PredictedProb)) +
geom_ribbon(aes(ymin = LL, ymax = UL, fill = Impacts), alpha = 0.2) +
geom_line(aes(colour = Impacts), size = 1)
# Predicted probabilities of migrants leaving sector specific, holding others at their means
contrasts(data$Sector_Company)
## MS SS
## AS 0 0
## MS 1 0
## SS 0 1
lr1 <- glm(formula = MW_left_covid ~ Total_Employees + Sector_Company + Age, data = data, family = "binomial")
newdata1 <- with(data, data.frame(Age = mean(Age), Total_Employees = mean(Total_Employees), Sector_Company = factor(c("AS", "MS", "SS"))))
newdata1$P <- predict(lr1, newdata = newdata1, type = "response")
newdata1
## Age Total_Employees Sector_Company P
## 1 14.6195 21.40671 AS 0.3777910
## 2 14.6195 21.40671 MS 0.3640189
## 3 14.6195 21.40671 SS 0.1651179
#Sector specific
newdata2 <- with(data, data.frame(Age = rep(seq(from = 6, to = 90, length.out = 15), 3),
Total_Employees = mean(Total_Employees), Sector_Company = factor(rep(c("AS", "MS", "SS"), each = 15))))
newdata2
## Age Total_Employees Sector_Company
## 1 6 21.40671 AS
## 2 12 21.40671 AS
## 3 18 21.40671 AS
## 4 24 21.40671 AS
## 5 30 21.40671 AS
## 6 36 21.40671 AS
## 7 42 21.40671 AS
## 8 48 21.40671 AS
## 9 54 21.40671 AS
## 10 60 21.40671 AS
## 11 66 21.40671 AS
## 12 72 21.40671 AS
## 13 78 21.40671 AS
## 14 84 21.40671 AS
## 15 90 21.40671 AS
## 16 6 21.40671 MS
## 17 12 21.40671 MS
## 18 18 21.40671 MS
## 19 24 21.40671 MS
## 20 30 21.40671 MS
## 21 36 21.40671 MS
## 22 42 21.40671 MS
## 23 48 21.40671 MS
## 24 54 21.40671 MS
## 25 60 21.40671 MS
## 26 66 21.40671 MS
## 27 72 21.40671 MS
## 28 78 21.40671 MS
## 29 84 21.40671 MS
## 30 90 21.40671 MS
## 31 6 21.40671 SS
## 32 12 21.40671 SS
## 33 18 21.40671 SS
## 34 24 21.40671 SS
## 35 30 21.40671 SS
## 36 36 21.40671 SS
## 37 42 21.40671 SS
## 38 48 21.40671 SS
## 39 54 21.40671 SS
## 40 60 21.40671 SS
## 41 66 21.40671 SS
## 42 72 21.40671 SS
## 43 78 21.40671 SS
## 44 84 21.40671 SS
## 45 90 21.40671 SS
newdata3 <- cbind(newdata2, predict(lr1, newdata = newdata2, type = "link", se = TRUE))
newdata3 <- within(newdata3, {
PredictedProb <- plogis(fit)
LL <- plogis(fit - (1.96 * se.fit))
UL <- plogis(fit + (1.96 * se.fit))
})
newdata3
## Age Total_Employees Sector_Company fit se.fit residual.scale
## 1 6 21.40671 AS -0.69044568 0.17061039 1
## 2 12 21.40671 AS -0.55713577 0.15843417 1
## 3 18 21.40671 AS -0.42382585 0.15368947 1
## 4 24 21.40671 AS -0.29051594 0.15705128 1
## 5 30 21.40671 AS -0.15720602 0.16803374 1
## 6 36 21.40671 AS -0.02389610 0.18528669 1
## 7 42 21.40671 AS 0.10941381 0.20725000 1
## 8 48 21.40671 AS 0.24272373 0.23259310 1
## 9 54 21.40671 AS 0.37603364 0.26033080 1
## 10 60 21.40671 AS 0.50934356 0.28977626 1
## 11 66 21.40671 AS 0.64265347 0.32045908 1
## 12 72 21.40671 AS 0.77596339 0.35205589 1
## 13 78 21.40671 AS 0.90927330 0.38434133 1
## 14 84 21.40671 AS 1.04258322 0.41715555 1
## 15 90 21.40671 AS 1.17589314 0.45038299 1
## 16 6 21.40671 MS -0.74947400 0.09795524 1
## 17 12 21.40671 MS -0.61616408 0.08483377 1
## 18 18 21.40671 MS -0.48285417 0.08557974 1
## 19 24 21.40671 MS -0.34954425 0.09988293 1
## 20 30 21.40671 MS -0.21623433 0.12310521 1
## 21 36 21.40671 MS -0.08292442 0.15119117 1
## 22 42 21.40671 MS 0.05038550 0.18190171 1
## 23 48 21.40671 MS 0.18369541 0.21411043 1
## 24 54 21.40671 MS 0.31700533 0.24723250 1
## 25 60 21.40671 MS 0.45031524 0.28094505 1
## 26 66 21.40671 MS 0.58362516 0.31505860 1
## 27 72 21.40671 MS 0.71693507 0.34945573 1
## 28 78 21.40671 MS 0.85024499 0.38406025 1
## 29 84 21.40671 MS 0.98355491 0.41882076 1
## 30 90 21.40671 MS 1.11686482 0.45370140 1
## 31 6 21.40671 SS -1.81214154 0.26235790 1
## 32 12 21.40671 SS -1.67883163 0.25840891 1
## 33 18 21.40671 SS -1.54552171 0.25931496 1
## 34 24 21.40671 SS -1.41221180 0.26502625 1
## 35 30 21.40671 SS -1.27890188 0.27524383 1
## 36 36 21.40671 SS -1.14559196 0.28949093 1
## 37 42 21.40671 SS -1.01228205 0.30720746 1
## 38 48 21.40671 SS -0.87897213 0.32783140 1
## 39 54 21.40671 SS -0.74566222 0.35085041 1
## 40 60 21.40671 SS -0.61235230 0.37582467 1
## 41 66 21.40671 SS -0.47904239 0.40239029 1
## 42 72 21.40671 SS -0.34573247 0.43025258 1
## 43 78 21.40671 SS -0.21242256 0.45917558 1
## 44 84 21.40671 SS -0.07911264 0.48897109 1
## 45 90 21.40671 SS 0.05419728 0.51948900 1
## UL LL PredictedProb
## 1 0.4119163 0.26408530 0.3339339
## 2 0.4386594 0.29574004 0.3642104
## 3 0.4693897 0.32628032 0.3956016
## 4 0.5043260 0.35472438 0.4278776
## 5 0.5429291 0.38070613 0.4607792
## 6 0.5840122 0.40442572 0.4940263
## 7 0.6261239 0.42634087 0.5273262
## 8 0.6678787 0.44691117 0.5603848
## 9 0.7081223 0.46649660 0.5929161
## 10 0.7459836 0.48534972 0.6246526
## 11 0.7808717 0.50363835 0.6553530
## 12 0.8124476 0.52147025 0.6848095
## 13 0.8405843 0.53891223 0.7128514
## 14 0.8653212 0.55600361 0.7393481
## 15 0.8868181 0.57276529 0.7642086
## 16 0.3641303 0.28060426 0.3209359
## 17 0.3893869 0.31379462 0.3506544
## 18 0.4218660 0.34285649 0.3815784
## 19 0.4616321 0.36695209 0.4134929
## 20 0.5062626 0.38757418 0.4461511
## 21 0.5531510 0.40630560 0.4792808
## 22 0.6003474 0.42405675 0.5125937
## 23 0.6464228 0.44128193 0.5457951
## 24 0.6903126 0.45820516 0.5785943
## 25 0.7312488 0.47493676 0.6107142
## 26 0.7687275 0.49152839 0.6419011
## 27 0.8024802 0.50799978 0.6719317
## 28 0.8324377 0.52435244 0.7006185
## 29 0.8586890 0.54057712 0.7278130
## 30 0.8814381 0.55665812 0.7534067
## 31 0.2145153 0.08896299 0.1403795
## 32 0.2364304 0.10107715 0.1572502
## 33 0.2616782 0.11367109 0.1757340
## 34 0.2905405 0.12656660 0.1958854
## 35 0.3231301 0.13962840 0.2177372
## 36 0.3593492 0.15277576 0.2412951
## 37 0.3988749 0.16597743 0.2665335
## 38 0.4411681 0.17923755 0.2933908
## 39 0.4855052 0.19258051 0.3217672
## 40 0.5310261 0.20603903 0.3515228
## 41 0.5767980 0.21964654 0.3824783
## 42 0.6218864 0.23343291 0.4144177
## 43 0.6654243 0.24742248 0.4470932
## 44 0.7066711 0.26163344 0.4802321
## 45 0.7450522 0.27607781 0.5135460
ggplot(newdata3, aes(x = Age, y = PredictedProb)) +
geom_ribbon(aes(ymin = LL, ymax = UL, fill = Sector_Company), alpha = 0.2) +
geom_line(aes(colour = Sector_Company), size = 1)
# Predicted probabilities of migrants leaving location specific, holding others at their means
contrasts(data$R_U)
## Urban
## Rural 0
## Urban 1
lr1 <- glm(formula = MW_left_covid ~ Total_Employees + R_U + Age, data = data, family = "binomial")
newdata1 <- with(data, data.frame(Age = mean(Age), Total_Employees = mean(Total_Employees), R_U = factor(c("Rural", "Urban"))))
newdata1$P <- predict(lr1, newdata = newdata1, type = "response")
newdata1
## Age Total_Employees R_U P
## 1 14.6195 21.40671 Rural 0.5418707
## 2 14.6195 21.40671 Urban 0.3318993
#Rural Vs Urban
newdata2 <- with(data, data.frame(Age = rep(seq(from = 6, to = 90, length.out = 15), 2),
Total_Employees = mean(Total_Employees), R_U = factor(rep(c("Rural", "Urban"), each = 15))))
newdata2
## Age Total_Employees R_U
## 1 6 21.40671 Rural
## 2 12 21.40671 Rural
## 3 18 21.40671 Rural
## 4 24 21.40671 Rural
## 5 30 21.40671 Rural
## 6 36 21.40671 Rural
## 7 42 21.40671 Rural
## 8 48 21.40671 Rural
## 9 54 21.40671 Rural
## 10 60 21.40671 Rural
## 11 66 21.40671 Rural
## 12 72 21.40671 Rural
## 13 78 21.40671 Rural
## 14 84 21.40671 Rural
## 15 90 21.40671 Rural
## 16 6 21.40671 Urban
## 17 12 21.40671 Urban
## 18 18 21.40671 Urban
## 19 24 21.40671 Urban
## 20 30 21.40671 Urban
## 21 36 21.40671 Urban
## 22 42 21.40671 Urban
## 23 48 21.40671 Urban
## 24 54 21.40671 Urban
## 25 60 21.40671 Urban
## 26 66 21.40671 Urban
## 27 72 21.40671 Urban
## 28 78 21.40671 Urban
## 29 84 21.40671 Urban
## 30 90 21.40671 Urban
newdata3 <- cbind(newdata2, predict(lr1, newdata = newdata2, type = "link", se = TRUE))
newdata3 <- within(newdata3, {
PredictedProb <- plogis(fit)
LL <- plogis(fit - (1.96 * se.fit))
UL <- plogis(fit + (1.96 * se.fit))
})
newdata3
## Age Total_Employees R_U fit se.fit residual.scale UL
## 1 6 21.40671 Rural -0.03769312 0.29451897 1 0.6317110
## 2 12 21.40671 Rural 0.10540284 0.29044086 1 0.6625474
## 3 18 21.40671 Rural 0.24849880 0.29056657 1 0.6938136
## 4 24 21.40671 Rural 0.39159476 0.29489073 1 0.7250359
## 5 30 21.40671 Rural 0.53469072 0.30323378 1 0.7556596
## 6 36 21.40671 Rural 0.67778668 0.31527683 1 0.7851153
## 7 42 21.40671 Rural 0.82088264 0.33061580 1 0.8128886
## 8 48 21.40671 Rural 0.96397859 0.34881614 1 0.8385743
## 9 54 21.40671 Rural 1.10707455 0.36945523 1 0.8619054
## 10 60 21.40671 Rural 1.25017051 0.39214819 1 0.8827549
## 11 66 21.40671 Rural 1.39326647 0.41655950 1 0.9011193
## 12 72 21.40671 Rural 1.53636243 0.44240478 1 0.9170920
## 13 78 21.40671 Rural 1.67945839 0.46944727 1 0.9308342
## 14 84 21.40671 Rural 1.82255435 0.49749177 1 0.9425481
## 15 90 21.40671 Rural 1.96565030 0.52637815 1 0.9524543
## 16 6 21.40671 Urban -0.90517656 0.09118390 1 0.3259729
## 17 12 21.40671 Urban -0.76208060 0.07524491 1 0.3510105
## 18 18 21.40671 Urban -0.61898464 0.07394363 1 0.3836569
## 19 24 21.40671 Urban -0.47588868 0.08793234 1 0.4246920
## 20 30 21.40671 Urban -0.33279272 0.11160243 1 0.4715179
## 21 36 21.40671 Urban -0.18969676 0.14013094 1 0.5212272
## 22 42 21.40671 Urban -0.04660081 0.17110473 1 0.5716936
## 23 48 21.40671 Urban 0.09649515 0.20340980 1 0.6213256
## 24 54 21.40671 Urban 0.23959111 0.23650124 1 0.6688822
## 25 60 21.40671 Urban 0.38268707 0.27009017 1 0.7134223
## 26 66 21.40671 Urban 0.52578303 0.30401173 1 0.7542939
## 27 72 21.40671 Urban 0.66887899 0.33816584 1 0.7911191
## 28 78 21.40671 Urban 0.81197495 0.37248854 1 0.8237629
## 29 84 21.40671 Urban 0.95507090 0.40693716 1 0.8522890
## 30 90 21.40671 Urban 1.09816686 0.44148224 1 0.8769082
## LL PredictedProb
## 1 0.3509308 0.4905778
## 2 0.3860702 0.5263263
## 3 0.4204292 0.5618070
## 4 0.4535367 0.5966665
## 5 0.4850925 0.6305765
## 6 0.5149566 0.6632445
## 7 0.5431116 0.6944237
## 8 0.5696195 0.7239177
## 9 0.5945826 0.7515833
## 10 0.6181162 0.7773294
## 11 0.6403330 0.8011132
## 12 0.6613350 0.8229353
## 13 0.6812108 0.8428328
## 14 0.7000362 0.8608723
## 15 0.7178758 0.8771431
## 16 0.2527693 0.2879879
## 17 0.2870898 0.3181947
## 18 0.3177971 0.3500124
## 19 0.3433871 0.3832234
## 20 0.3655087 0.4175613
## 21 0.3859536 0.4527175
## 22 0.4056528 0.4883519
## 23 0.4250227 0.5241051
## 24 0.4442450 0.5596129
## 25 0.4633932 0.5945210
## 26 0.4824872 0.6284990
## 27 0.5015185 0.6612521
## 28 0.5204629 0.6925302
## 29 0.5392874 0.7221338
## 30 0.5579538 0.7499165
ggplot(newdata3, aes(x = Age, y = PredictedProb)) +
geom_ribbon(aes(ymin = LL, ymax = UL, fill = R_U), alpha = 0.2) +
geom_line(aes(colour = R_U), size = 1)
#Ordinal Logistic Regression (ordered)
######################################################
# relevel - only for unordered factors --
library(foreign)
library(MASS)
## Warning: package 'MASS' was built under R version 4.0.3
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
library(Hmisc)
## Warning: package 'Hmisc' was built under R version 4.0.3
## Loading required package: lattice
## Loading required package: survival
##
## Attaching package: 'survival'
## The following object is masked from 'package:aod':
##
## rats
## Loading required package: Formula
## Warning: package 'Formula' was built under R version 4.0.3
##
## Attaching package: 'Hmisc'
## The following object is masked from 'package:sjmisc':
##
## %nin%
## The following objects are masked from 'package:dplyr':
##
## src, summarize
## The following objects are masked from 'package:plyr':
##
## is.discrete, summarize
## The following objects are masked from 'package:base':
##
## format.pval, units
library(reshape2)
##
## Attaching package: 'reshape2'
## The following object is masked from 'package:tidyr':
##
## smiths
ggplot(data, aes(x = Impacts, y = prc_loss)) +
coord_cartesian(ylim = c(0, 100)) +
geom_boxplot(size = .75) +
facet_grid(R_U ~ Sector_Company, margins = TRUE) +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1))
## Warning: Removed 4 rows containing non-finite values (stat_boxplot).
levels(data$impact_severity)
## [1] "None" "Marginal" "Moderate" "Severe"
olr <- polr(formula = impact_severity ~ prc_loss + Age + R_U + Sector_Company + Size_Company + Type_Company + Total_Employees + Migrant_Employees, data = data, Hess = TRUE)
## Warning in polr(formula = impact_severity ~ prc_loss + Age + R_U +
## Sector_Company + : design appears to be rank-deficient, so dropping some coefs
summary(olr)
## Call:
## polr(formula = impact_severity ~ prc_loss + Age + R_U + Sector_Company +
## Size_Company + Type_Company + Total_Employees + Migrant_Employees,
## data = data, Hess = TRUE)
##
## Coefficients:
## Value Std. Error t value
## prc_loss 0.035996 0.003892 9.2495
## Age 0.009306 0.006550 1.4209
## R_UUrban -0.526278 0.334096 -1.5752
## Sector_CompanyMS 3.549754 0.289653 12.2552
## Sector_CompanySS 3.438334 0.348203 9.8745
## Size_Company.L -0.146982 1.070047 -0.1374
## Size_Company.Q -2.799971 1.361495 -2.0565
## Type_CompanyCLS_Pvt 2.345108 0.720645 3.2542
## Type_CompanyCLS_Pub -0.576508 1.169628 -0.4929
## Type_CompanyUL 0.627570 0.710811 0.8829
## Total_Employees -0.007613 0.002550 -2.9853
## Migrant_Employees 0.022307 0.007498 2.9749
##
## Intercepts:
## Value Std. Error t value
## None|Marginal 2.3290 1.0464 2.2258
## Marginal|Moderate 7.0735 1.0741 6.5855
## Moderate|Severe 11.4404 1.1286 10.1371
##
## Residual Deviance: 1294.512
## AIC: 1324.512
## (35 observations deleted due to missingness)
ctable <- coef(summary(olr))
tab_model(olr)
|
|
impact_severity
|
|
Predictors
|
Odds Ratios
|
CI
|
p
|
|
None|Marginal
|
10.27
|
10.19 – 10.35
|
0.026
|
|
Marginal|Moderate
|
1180.29
|
1165.22 – 1195.56
|
<0.001
|
|
Moderate|Severe
|
93003.53
|
48276.43 – 179169.32
|
<0.001
|
|
prc_loss
|
1.04
|
0.59 – 1.83
|
<0.001
|
|
Age
|
1.01
|
0.51 – 2.00
|
0.156
|
|
R_U [Urban]
|
0.59
|
0.07 – 4.82
|
0.116
|
|
Sector_Company [MS]
|
34.80
|
2.41 – 503.62
|
<0.001
|
|
Sector_Company [SS]
|
31.14
|
7.57 – 128.08
|
<0.001
|
|
Size_Company [linear]
|
0.86
|
0.09 – 8.57
|
0.891
|
|
Size_Company [quadratic]
|
0.06
|
0.02 – 0.25
|
0.040
|
|
Type_Company [CLS_Pvt]
|
10.43
|
10.38 – 10.49
|
0.001
|
|
Type_Company [CLS_Pub]
|
0.56
|
0.55 – 0.57
|
0.622
|
|
Type_Company [UL]
|
1.87
|
0.24 – 14.60
|
0.378
|
|
Total_Employees
|
0.99
|
0.12 – 8.17
|
0.003
|
|
Migrant_Employees
|
1.02
|
0.11 – 9.37
|
0.003
|
|
Observations
|
919
|
|
R2 Nagelkerke
|
0.450
|
# add p value
p <- pnorm(abs(ctable[, "t value"]), lower.tail = FALSE) * 2
ctable <- cbind(ctable, "p value" = p)
View(ctable)
# add confidence interval
#ci <- confint(olr)
#ci <- confint.default(olr)
#ci
#exp(coef(olr))
#exp(cbind(OR = coef(olr), ci))
#predicted probabilities
olr1 <- polr(formula = impact_severity ~ MW_left_covid + prc_loss + Age + Sector_Company + Total_Employees + Migrant_Employees, data = data, Hess = TRUE)
summary(olr1)
## Call:
## polr(formula = impact_severity ~ MW_left_covid + prc_loss + Age +
## Sector_Company + Total_Employees + Migrant_Employees, data = data,
## Hess = TRUE)
##
## Coefficients:
## Value Std. Error t value
## MW_left_covidYes 0.380439 0.153061 2.486
## prc_loss 0.034253 0.003588 9.547
## Age 0.010056 0.006027 1.668
## Sector_CompanyMS 1.840401 0.221361 8.314
## Sector_CompanySS 2.134606 0.286398 7.453
## Total_Employees -0.002403 0.002191 -1.097
## Migrant_Employees 0.007497 0.006340 1.182
##
## Intercepts:
## Value Std. Error t value
## None|Marginal -0.6360 0.2491 -2.5536
## Marginal|Moderate 3.7603 0.2880 13.0587
## Moderate|Severe 7.6643 0.3804 20.1494
##
## Residual Deviance: 1474.894
## AIC: 1494.894
## (12 observations deleted due to missingness)