setwd("C:/Users/bha/Desktop/caroline/Data")
library(nnet)
## Warning: package 'nnet' was built under R version 4.0.5
library(MASS)
load('dat2.RData')
### Logit Model
m1_logit <- glm(bin ~ AGE + factor(Gend3) + factor(Disab3) +
factor(IMD10) + factor(BMIG) + factor(WorkStat10) +
factor(UrbRur2) + factor(CoastComm) + factor(Educ6)+factor(NSSEC5), data = dat2, family = 'binomial')
summary(m1_logit)
##
## Call:
## glm(formula = bin ~ AGE + factor(Gend3) + factor(Disab3) + factor(IMD10) +
## factor(BMIG) + factor(WorkStat10) + factor(UrbRur2) + factor(CoastComm) +
## factor(Educ6) + factor(NSSEC5), family = "binomial", data = dat2)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.2075 -1.1338 0.6640 0.9448 1.7805
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 2.7351644 1.0178921 2.687 0.00721 **
## AGE -0.0008142 0.0067722 -0.120 0.90431
## factor(Gend3)2 -0.3406836 0.1730423 -1.969 0.04898 *
## factor(Disab3)2.00 0.3919080 0.2714182 1.444 0.14876
## factor(Disab3)3.00 0.5689855 0.2374167 2.397 0.01655 *
## factor(IMD10)2 -0.8531074 0.5071475 -1.682 0.09254 .
## factor(IMD10)3 -1.1764853 0.5199868 -2.263 0.02366 *
## factor(IMD10)4 -0.8021813 0.5047724 -1.589 0.11202
## factor(IMD10)5 -1.0550014 0.5271769 -2.001 0.04537 *
## factor(IMD10)6 -0.6120645 0.5068226 -1.208 0.22718
## factor(IMD10)7 -1.0935183 0.5033553 -2.172 0.02982 *
## factor(IMD10)8 -1.1839899 0.4616395 -2.565 0.01033 *
## factor(IMD10)9 -1.0170773 0.4708700 -2.160 0.03077 *
## factor(IMD10)10 -1.2325669 0.4796267 -2.570 0.01017 *
## factor(BMIG)2 0.0245081 0.6943862 0.035 0.97184
## factor(BMIG)3 -0.2852124 0.6960783 -0.410 0.68200
## factor(BMIG)4 -0.7456920 0.7077038 -1.054 0.29203
## factor(BMIG)5 -1.1308549 0.8630103 -1.310 0.19007
## factor(WorkStat10)2 -0.4795079 0.2561335 -1.872 0.06119 .
## factor(WorkStat10)5 -0.5552338 0.2151930 -2.580 0.00988 **
## factor(UrbRur2)2 0.2209164 0.1739597 1.270 0.20411
## factor(CoastComm)2 -0.6446887 0.4381490 -1.471 0.14119
## factor(Educ6)2 -0.2395086 0.2364180 -1.013 0.31103
## factor(Educ6)3 -0.4920173 0.2104821 -2.338 0.01941 *
## factor(NSSEC5)2.00 -0.5337077 0.2219062 -2.405 0.01617 *
## factor(NSSEC5)3.00 0.0304553 0.2541382 0.120 0.90461
## factor(NSSEC5)4.00 0.1241942 0.3811848 0.326 0.74457
## factor(NSSEC5)5.00 -0.5352807 0.3665898 -1.460 0.14425
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 984.36 on 735 degrees of freedom
## Residual deviance: 893.25 on 708 degrees of freedom
## AIC: 949.25
##
## Number of Fisher Scoring iterations: 4
## odds ratios only
exp(coef(m1_logit))
## (Intercept) AGE factor(Gend3)2 factor(Disab3)2.00
## 15.4122773 0.9991862 0.7112839 1.4798015
## factor(Disab3)3.00 factor(IMD10)2 factor(IMD10)3 factor(IMD10)4
## 1.7664740 0.4260889 0.3083606 0.4483499
## factor(IMD10)5 factor(IMD10)6 factor(IMD10)7 factor(IMD10)8
## 0.3481920 0.5422303 0.3350357 0.3060552
## factor(IMD10)9 factor(IMD10)10 factor(BMIG)2 factor(BMIG)3
## 0.3616504 0.2915433 1.0248109 0.7518545
## factor(BMIG)4 factor(BMIG)5 factor(WorkStat10)2 factor(WorkStat10)5
## 0.4744059 0.3227572 0.6190880 0.5739380
## factor(UrbRur2)2 factor(CoastComm)2 factor(Educ6)2 factor(Educ6)3
## 1.2472191 0.5248259 0.7870145 0.6113918
## factor(NSSEC5)2.00 factor(NSSEC5)3.00 factor(NSSEC5)4.00 factor(NSSEC5)5.00
## 0.5864266 1.0309238 1.1322357 0.5855049
## odds ratios and 95% CI
exp(cbind(OR = coef(m1_logit), confint(m1_logit)))
## Waiting for profiling to be done...
## OR 2.5 % 97.5 %
## (Intercept) 15.4122773 2.18995112 121.4182195
## AGE 0.9991862 0.98598464 1.0125532
## factor(Gend3)2 0.7112839 0.50591130 0.9975482
## factor(Disab3)2.00 1.4798015 0.87049719 2.5261972
## factor(Disab3)3.00 1.7664740 1.10977284 2.8188602
## factor(IMD10)2 0.4260889 0.15088961 1.1211035
## factor(IMD10)3 0.3083606 0.10652933 0.8313478
## factor(IMD10)4 0.4483499 0.15922171 1.1711386
## factor(IMD10)5 0.3481920 0.11894188 0.9554762
## factor(IMD10)6 0.5422303 0.19220601 1.4258505
## factor(IMD10)7 0.3350357 0.11909255 0.8704558
## factor(IMD10)8 0.3060552 0.11691653 0.7268869
## factor(IMD10)9 0.3616504 0.13610327 0.8771814
## factor(IMD10)10 0.2915433 0.10798729 0.7197175
## factor(BMIG)2 1.0248109 0.23788948 3.8545560
## factor(BMIG)3 0.7518545 0.17406723 2.8387089
## factor(BMIG)4 0.4744059 0.10776693 1.8344964
## factor(BMIG)5 0.3227572 0.05505913 1.6893888
## factor(WorkStat10)2 0.6190880 0.37470207 1.0242207
## factor(WorkStat10)5 0.5739380 0.37539918 0.8735443
## factor(UrbRur2)2 1.2472191 0.88724178 1.7556540
## factor(CoastComm)2 0.5248259 0.21555625 1.2161073
## factor(Educ6)2 0.7870145 0.49599941 1.2547809
## factor(Educ6)3 0.6113918 0.40410207 0.9230238
## factor(NSSEC5)2.00 0.5864266 0.37928226 0.9061259
## factor(NSSEC5)3.00 1.0309238 0.62814888 1.7033084
## factor(NSSEC5)4.00 1.1322357 0.54230449 2.4342753
## factor(NSSEC5)5.00 0.5855049 0.28418081 1.2001136
### Multinomial baseline logit model
dat2$tri <- as.factor(dat2$tri)
dat2$tri2 <- relevel(dat2$tri, ref = 3)
library(nnet)
m1_multinom <- multinom(tri2 ~ AGE + factor(Gend3) + factor(Disab3) +
factor(IMD10) + factor(BMIG) + factor(WorkStat10) +
factor(UrbRur2) + factor(CoastComm) + factor(Educ6)+factor(NSSEC5), data = dat2)
## # weights: 87 (56 variable)
## initial value 808.578644
## iter 10 value 634.283470
## iter 20 value 609.165895
## iter 30 value 607.248752
## iter 40 value 606.390480
## iter 50 value 606.301745
## iter 60 value 606.227498
## iter 70 value 606.223227
## final value 606.223173
## converged
summary(m1_multinom)
## Call:
## multinom(formula = tri2 ~ AGE + factor(Gend3) + factor(Disab3) +
## factor(IMD10) + factor(BMIG) + factor(WorkStat10) + factor(UrbRur2) +
## factor(CoastComm) + factor(Educ6) + factor(NSSEC5), data = dat2)
##
## Coefficients:
## (Intercept) AGE factor(Gend3)2 factor(Disab3)2.00 factor(Disab3)3.00
## 1 2.145169 -0.008469122 -0.1590435 0.6394732 0.7351749
## 2 -16.799887 -0.019905601 0.4744624 0.8200478 0.5639266
## factor(IMD10)2 factor(IMD10)3 factor(IMD10)4 factor(IMD10)5 factor(IMD10)6
## 1 -0.4962602 -0.7878149 -0.2465553 -0.2299663 0.09042603
## 2 1.4577654 1.5683910 1.9836527 2.5912171 2.32673095
## factor(IMD10)7 factor(IMD10)8 factor(IMD10)9 factor(IMD10)10 factor(BMIG)2
## 1 -0.1597315 -0.6917853 -0.530486 -0.9634356 0.6349464
## 2 2.8002388 1.8281893 1.773263 0.9363959 13.5801667
## factor(BMIG)3 factor(BMIG)4 factor(BMIG)5 factor(WorkStat10)2
## 1 0.1981778 -0.08704732 -0.3275833 -0.487005757
## 2 13.2601529 13.73637318 14.0335956 -0.007647693
## factor(WorkStat10)5 factor(UrbRur2)2 factor(CoastComm)2 factor(Educ6)2
## 1 -0.4988571 0.04793871 -0.2304298 -0.4727966
## 2 0.1518663 -0.48793159 1.8340285 -0.5483517
## factor(Educ6)3 factor(NSSEC5)2.00 factor(NSSEC5)3.00 factor(NSSEC5)4.00
## 1 -0.7080472 -0.4565365 -0.1812350 0.04880388
## 2 -0.5890083 0.2031513 -0.6359764 -0.18404419
## factor(NSSEC5)5.00
## 1 -0.5791999
## 2 -0.2197591
##
## Std. Errors:
## (Intercept) AGE factor(Gend3)2 factor(Disab3)2.00 factor(Disab3)3.00
## 1 1.087917 0.008066969 0.1995387 0.3043036 0.2622030
## 2 1.362008 0.010835389 0.2833641 0.4231386 0.3751397
## factor(IMD10)2 factor(IMD10)3 factor(IMD10)4 factor(IMD10)5 factor(IMD10)6
## 1 0.5514726 0.5674838 0.5573031 0.6104916 0.5725756
## 2 1.1939055 1.1844059 1.1809883 1.2034294 1.1799432
## factor(IMD10)7 factor(IMD10)8 factor(IMD10)9 factor(IMD10)10 factor(BMIG)2
## 1 0.5738518 0.4984278 0.5113751 0.5127815 0.7098938
## 2 1.1696752 1.1295325 1.1371590 1.1662045 0.4262360
## factor(BMIG)3 factor(BMIG)4 factor(BMIG)5 factor(WorkStat10)2
## 1 0.7096654 0.7255603 0.9276842 0.3072961
## 2 0.4280840 0.4406446 0.6302095 0.4075177
## factor(WorkStat10)5 factor(UrbRur2)2 factor(CoastComm)2 factor(Educ6)2
## 1 0.2534737 0.2037399 0.4731469 0.2838423
## 2 0.3470588 0.2750452 1.1037827 0.3766972
## factor(Educ6)3 factor(NSSEC5)2.00 factor(NSSEC5)3.00 factor(NSSEC5)4.00
## 1 0.2450305 0.2619431 0.2891527 0.4439384
## 2 0.3395679 0.3418254 0.4316972 0.6389025
## factor(NSSEC5)5.00
## 1 0.4040604
## 2 0.5957776
##
## Residual Deviance: 1212.446
## AIC: 1324.446
z <- summary(m1_multinom)$coefficients/summary(m1_multinom)$standard.errors
z
## (Intercept) AGE factor(Gend3)2 factor(Disab3)2.00 factor(Disab3)3.00
## 1 1.971812 -1.049852 -0.7970559 2.101431 2.803838
## 2 -12.334646 -1.837092 1.6743915 1.938013 1.503244
## factor(IMD10)2 factor(IMD10)3 factor(IMD10)4 factor(IMD10)5 factor(IMD10)6
## 1 -0.899882 -1.388260 -0.4424079 -0.3766903 0.1579286
## 2 1.221006 1.324201 1.6796549 2.1531942 1.9719009
## factor(IMD10)7 factor(IMD10)8 factor(IMD10)9 factor(IMD10)10 factor(BMIG)2
## 1 -0.2783497 -1.387935 -1.037372 -1.8788422 0.8944244
## 2 2.3940312 1.618536 1.559380 0.8029431 31.8606775
## factor(BMIG)3 factor(BMIG)4 factor(BMIG)5 factor(WorkStat10)2
## 1 0.2792553 -0.1199725 -0.3531194 -1.58480933
## 2 30.9755863 31.1733638 22.2681451 -0.01876653
## factor(WorkStat10)5 factor(UrbRur2)2 factor(CoastComm)2 factor(Educ6)2
## 1 -1.9680824 0.2352936 -0.4870153 -1.665701
## 2 0.4375809 -1.7740053 1.6615848 -1.455683
## factor(Educ6)3 factor(NSSEC5)2.00 factor(NSSEC5)3.00 factor(NSSEC5)4.00
## 1 -2.889629 -1.7428841 -0.6267796 0.1099339
## 2 -1.734582 0.5943131 -1.4732003 -0.2880630
## factor(NSSEC5)5.00
## 1 -1.4334489
## 2 -0.3688609
# 2-tailed z test
p <- (1 - pnorm(abs(z), 0, 1)) * 2
p
## (Intercept) AGE factor(Gend3)2 factor(Disab3)2.00 factor(Disab3)3.00
## 1 0.04863105 0.29378629 0.42541857 0.03560311 0.005049825
## 2 0.00000000 0.06619638 0.09405367 0.05262169 0.132776028
## factor(IMD10)2 factor(IMD10)3 factor(IMD10)4 factor(IMD10)5 factor(IMD10)6
## 1 0.3681831 0.1650580 0.65819405 0.70640376 0.87451310
## 2 0.2220839 0.1854364 0.09302449 0.03130342 0.04862092
## factor(IMD10)7 factor(IMD10)8 factor(IMD10)9 factor(IMD10)10 factor(BMIG)2
## 1 0.78074394 0.1651569 0.2995627 0.06026604 0.3710949
## 2 0.01666433 0.1055471 0.1189064 0.42200760 0.0000000
## factor(BMIG)3 factor(BMIG)4 factor(BMIG)5 factor(WorkStat10)2
## 1 0.7800489 0.9045049 0.7239989 0.1130097
## 2 0.0000000 0.0000000 0.0000000 0.9850274
## factor(WorkStat10)5 factor(UrbRur2)2 factor(CoastComm)2 factor(Educ6)2
## 1 0.04905856 0.81398084 0.62624751 0.09577291
## 2 0.66169010 0.07606227 0.09659605 0.14548031
## factor(Educ6)3 factor(NSSEC5)2.00 factor(NSSEC5)3.00 factor(NSSEC5)4.00
## 1 0.003856967 0.08135386 0.5308037 0.9124618
## 2 0.082814895 0.55230272 0.1406970 0.7732985
## factor(NSSEC5)5.00
## 1 0.1517296
## 2 0.7122314
exp(coef(m1_multinom)) ## odd ratio
## (Intercept) AGE factor(Gend3)2 factor(Disab3)2.00 factor(Disab3)3.00
## 1 8.543483e+00 0.9915666 0.8529592 1.895482 2.085847
## 2 5.057102e-08 0.9802912 1.6071499 2.270608 1.757560
## factor(IMD10)2 factor(IMD10)3 factor(IMD10)4 factor(IMD10)5 factor(IMD10)6
## 1 0.6088032 0.4548376 0.7814881 0.7945604 1.094641
## 2 4.2963483 4.7989207 7.2692467 13.3460047 10.244397
## factor(IMD10)7 factor(IMD10)8 factor(IMD10)9 factor(IMD10)10 factor(BMIG)2
## 1 0.8523726 0.5006814 0.588319 0.3815797 1.886921e+00
## 2 16.4485745 6.2226094 5.890043 2.5507715 7.902991e+05
## factor(BMIG)3 factor(BMIG)4 factor(BMIG)5 factor(WorkStat10)2
## 1 1.219179e+00 9.166337e-01 7.206633e-01 0.6144635
## 2 5.738670e+05 9.239130e+05 1.243693e+06 0.9923815
## factor(WorkStat10)5 factor(UrbRur2)2 factor(CoastComm)2 factor(Educ6)2
## 1 0.6072243 1.0491064 0.7941922 0.6232568
## 2 1.1640046 0.6138949 6.2590505 0.5779016
## factor(Educ6)3 factor(NSSEC5)2.00 factor(NSSEC5)3.00 factor(NSSEC5)4.00
## 1 0.4926052 0.6334739 0.8342393 1.050014
## 2 0.5548773 1.2252579 0.5294183 0.831899
## factor(NSSEC5)5.00
## 1 0.5603465
## 2 0.8027122
exp(confint(m1_multinom)) ## CI
## , , 1
##
## 2.5 % 97.5 %
## (Intercept) 1.0129733 72.0562857
## AGE 0.9760123 1.0073689
## factor(Gend3)2 0.5768735 1.2611767
## factor(Disab3)2.00 1.0439891 3.4414651
## factor(Disab3)3.00 1.2476557 3.4871452
## factor(IMD10)2 0.2065678 1.7942842
## factor(IMD10)3 0.1495592 1.3832460
## factor(IMD10)4 0.2621471 2.3296982
## factor(IMD10)5 0.2401461 2.6289255
## factor(IMD10)6 0.3563644 3.3623949
## factor(IMD10)7 0.2767999 2.6247810
## factor(IMD10)8 0.1884946 1.3299153
## factor(IMD10)9 0.2159382 1.6028620
## factor(IMD10)10 0.1396705 1.0424750
## factor(BMIG)2 0.4693437 7.5860629
## factor(BMIG)3 0.3033886 4.8993199
## factor(BMIG)4 0.2211046 3.8000902
## factor(BMIG)5 0.1169731 4.4399558
## factor(WorkStat10)2 0.3364536 1.1221916
## factor(WorkStat10)5 0.3694808 0.9979443
## factor(UrbRur2)2 0.7037133 1.5640235
## factor(CoastComm)2 0.3141827 2.0075622
## factor(Educ6)2 0.3573222 1.0871115
## factor(Educ6)3 0.3047395 0.7962865
## factor(NSSEC5)2.00 0.3791075 1.0585103
## factor(NSSEC5)3.00 0.4733293 1.4703403
## factor(NSSEC5)4.00 0.4398598 2.5065494
## factor(NSSEC5)5.00 0.2538165 1.2370678
##
## , , 2
##
## 2.5 % 97.5 %
## (Intercept) 3.503954e-09 7.298692e-07
## AGE 9.596923e-01 1.001332e+00
## factor(Gend3)2 9.222666e-01 2.800634e+00
## factor(Disab3)2.00 9.907545e-01 5.203774e+00
## factor(Disab3)3.00 8.425404e-01 3.666314e+00
## factor(IMD10)2 4.138522e-01 4.460193e+01
## factor(IMD10)3 4.709507e-01 4.890032e+01
## factor(IMD10)4 7.181752e-01 7.357808e+01
## factor(IMD10)5 1.261799e+00 1.411602e+02
## factor(IMD10)6 1.014184e+00 1.034799e+02
## factor(IMD10)7 1.661495e+00 1.628387e+02
## factor(IMD10)8 6.800060e-01 5.694195e+01
## factor(IMD10)9 6.341134e-01 5.471041e+01
## factor(IMD10)10 2.594159e-01 2.508110e+01
## factor(BMIG)2 3.427510e+05 1.822234e+06
## factor(BMIG)3 2.479850e+05 1.327997e+06
## factor(BMIG)4 3.895415e+05 2.191333e+06
## factor(BMIG)5 3.616412e+05 4.277090e+06
## factor(WorkStat10)2 4.464769e-01 2.205760e+00
## factor(WorkStat10)5 5.895732e-01 2.298115e+00
## factor(UrbRur2)2 3.580759e-01 1.052478e+00
## factor(CoastComm)2 7.193942e-01 5.445653e+01
## factor(Educ6)2 2.761904e-01 1.209203e+00
## factor(Educ6)3 2.852044e-01 1.079537e+00
## factor(NSSEC5)2.00 6.269966e-01 2.394362e+00
## factor(NSSEC5)3.00 2.271630e-01 1.233844e+00
## factor(NSSEC5)4.00 2.378131e-01 2.910083e+00
## factor(NSSEC5)5.00 2.497083e-01 2.580398e+00