require("knitr")
## Loading required package: knitr
library(car)
## Loading required package: carData
library(VGAM)
## Loading required package: stats4
## Loading required package: splines
## 
## Attaching package: 'VGAM'
## The following object is masked from 'package:car':
## 
##     logit
library(stargazer)
## 
## Please cite as:
##  Hlavac, Marek (2018). stargazer: Well-Formatted Regression and Summary Statistics Tables.
##  R package version 5.2.2. https://CRAN.R-project.org/package=stargazer
library(survey)
## Loading required package: grid
## Loading required package: Matrix
## Loading required package: survival
## 
## Attaching package: 'survey'
## The following object is masked from 'package:VGAM':
## 
##     calibrate
## The following object is masked from 'package:graphics':
## 
##     dotchart
library(ggplot2)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following object is masked from 'package:car':
## 
##     recode
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(haven)
NSDUH_2019 <- read_sav("NSDUH_2019.SAV")
View(NSDUH_2019)

1. Define an ordinal or multinomial outcome variable of your choosing and define how you will recode the original variable.

The variable DSTHOP30 will be used to how often a respondent felt hopeless within the last 30 days. The variable will have to be releveled in the following ways with 1:2 equalling 3 or feeling hopeless most of the time, 3:4=2 feeling hopeless some of the time, and 5=1 as not feeling hopeless at all. The item was essentially reverse recoded, as lower numbers were indicative of more hopelessness, this is now the opposite.

2. State a research question about what factors you believe will affect your outcome variable.

How might socio demographic variables such as race, education, marital status among others effect the likelihood of being depressed within the last 30 days? DO things such as being of low socioeconomic status, or of socially marginialized groups increase the likelihood of hopelessness occuring within the last 30 days.

## Main variable hopeless last thirty days
NSDUH_2019$hopelesslst30days<-Recode(NSDUH_2019$dsthop30,
                               recodes="1:2=3;3:4=2;5=1; else=NA", as.factor = T)
NSDUH_2019$hopelesslst30days<-relevel(NSDUH_2019$hopelesslst30days, ref = "1")
NSDUH_2019$hopelesslst30daysnum<-car::Recode(NSDUH_2019$dsthop30,
                                recodes="1:2=3;3:4=2;5=1; else=NA", as.factor = F)

## marital status
NSDUH_2019$marst<-Recode(NSDUH_2019$irmarit, recodes="1='married'; 2='divorced'; 3='widowed'; 4='separated'; else=NA", as.factor=T)
NSDUH_2019$marst<-relevel(NSDUH_2019$marst, ref='married')

## education recodes
NSDUH_2019$educ<-Recode(NSDUH_2019$ireduhighst2, recodes="1:7='LssThnHgh'; 8='highschool'; 9='someCollege'; 10='associates'; 11='colgrad';else=NA", as.factor=T)
NSDUH_2019$educ<-relevel(NSDUH_2019$educ, ref='colgrad')

## sexuality recodes
NSDUH_2019$sexuality<-Recode(NSDUH_2019$sexident, recodes="1='Heterosexual'; 2='Les/Gay'; 3='Bisexual';else=NA", as.factor=T)
NSDUH_2019$sexuality<-relevel(NSDUH_2019$sexuality, ref='Heterosexual')

## gender recodes
NSDUH_2019$male<-as.factor(ifelse(NSDUH_2019$irsex==1, "Male", "Female"))

## Race recoded items
NSDUH_2019$black<-Recode(NSDUH_2019$newrace2, recodes="2=1; 9=NA; else=0")
NSDUH_2019$white<-Recode(NSDUH_2019$newrace2, recodes="1=1; 9=NA; else=0")
NSDUH_2019$other<-Recode(NSDUH_2019$newrace2, recodes="3:4=1; 9=NA; else=0")
NSDUH_2019$mult_race<-Recode(NSDUH_2019$newrace2, recodes="6=1; 9=NA; else=0")
NSDUH_2019$asian<-Recode(NSDUH_2019$newrace2, recodes="5=1; 9=NA; else=0")
NSDUH_2019$hispanic<-Recode(NSDUH_2019$newrace2, recodes="7=1; 9=NA; else=0")
NSDUH_2019$race_eth<-Recode(NSDUH_2019$newrace2,
                          recodes="1='white'; 2='black'; 3='other'; 4='asian'; 5='mult_race'; 6='hispanic'; else=NA",
                          as.factor = T)
NSDUH_2019$race_eth<-relevel(NSDUH_2019$race_eth, ref='white')
NSDUH_2019$lst_alc_use2<-Recode(NSDUH_2019$iralcrc, recodes="1='last 30days'; 2='12>1month'; 3='>12months'; else=NA", as.factor=T)
NSDUH_2019$dep_year2<-Recode(NSDUH_2019$amdeyr, recodes="1=1; 2=0;else=NA", as.factor=T)
NSDUH_2019$age_cat<-Recode(NSDUH_2019$age2, recodes="7:8='18-19'; 9:10='20-21'; 11='22-23'; 12='24-25'; 13='26-29'; 14='30-34'; 15='35-49'; 16='50-64'; 17='65+'; else=NA", as.factor=T)
options(survey.lonely.psu = "adjust")

Ordinal Regression Model Work

sub<-NSDUH_2019%>%
  select(hopelesslst30days,hopelesslst30daysnum, sexuality, male,
         age_cat,race_eth, marst, educ,white, black, hispanic, asian, other, mult_race,
         lst_alc_use2, dep_year2, vestr,analwtc ) %>%
  filter( complete.cases(.))
  
options(survey.lonely.psu = "adjust")
des<-svydesign(ids=~1, strata=~vestr, weights=~analwtc, data =sub )
fit.solr1<-svyolr(hopelesslst30days~race_eth+educ+age_cat+sexuality+male+marst+lst_alc_use2+dep_year2,des)
summary(fit.solr1)
## Call:
## svyolr(hopelesslst30days ~ race_eth + educ + age_cat + sexuality + 
##     male + marst + lst_alc_use2 + dep_year2, des)
## 
## Coefficients:
##                               Value Std. Error     t value
## race_ethasian           -0.03929821 0.26220764  -0.1498744
## race_ethblack           -0.18913942 0.05250227  -3.6025003
## race_ethhispanic         0.20156861 0.10441323   1.9304892
## race_ethmult_race        0.17656969 0.08376240   2.1079826
## race_ethother           -0.05227946 0.19614616  -0.2665332
## educassociates           0.15808757 0.06672784   2.3691396
## educhighschool           0.24397346 0.05363801   4.5485184
## educLssThnHgh            0.43092803 0.07665501   5.6216554
## educsomeCollege          0.15243853 0.05099274   2.9894162
## age_cat20-21            -0.06813294 0.08426914  -0.8085159
## age_cat22-23            -0.05209686 0.08185702  -0.6364373
## age_cat24-25            -0.12491970 0.08151124  -1.5325457
## age_cat26-29            -0.18983909 0.08209346  -2.3124751
## age_cat30-34            -0.33104437 0.08423084  -3.9302040
## age_cat35-49            -0.54207064 0.08099116  -6.6929606
## age_cat50-64            -0.99101986 0.09110365 -10.8779376
## age_cat65+              -1.14567598 0.10071579 -11.3753367
## sexualityBisexual        0.48774346 0.06910503   7.0580025
## sexualityLes/Gay         0.20263788 0.12967809   1.5626224
## maleMale                -0.19199793 0.03806454  -5.0440097
## marstdivorced            0.23318698 0.11080016   2.1045727
## marstseparated           0.50798401 0.05114915   9.9314264
## marstwidowed             0.44667097 0.05942071   7.5170927
## lst_alc_use212>1month   -0.06237661 0.06583135  -0.9475214
## lst_alc_use2last 30days -0.16678440 0.05428399  -3.0724415
## dep_year21               2.28637049 0.05989444  38.1733373
## 
## Intercepts:
##     Value    Std. Error t value 
## 1|2   0.5957   0.1041     5.7233
## 2|3   3.4651   0.1104    31.3904
fit.solr1$deviance+2*length(fit.solr1$coefficients)
## [1] 37380.61
ex1<-svyglm(I(hopelesslst30daysnum>1)~race_eth+educ+age_cat+sexuality+male+marst+lst_alc_use2+dep_year2,des, family="binomial")
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!
ex2<-svyglm(I(hopelesslst30daysnum>2)~race_eth+educ+age_cat+sexuality+male+marst+lst_alc_use2+dep_year2,des, family="binomial")
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!
round(exp(rbind(coef(ex1)[-1],
                coef(ex2)[-1])),3)
##      race_ethasian race_ethblack race_ethhispanic race_ethmult_race
## [1,]         0.924         0.807            1.291             1.187
## [2,]         1.505         1.102            1.098             1.357
##      race_ethother educassociates educhighschool educLssThnHgh educsomeCollege
## [1,]         0.885          1.171          1.192         1.440           1.137
## [2,]         1.653          1.493          3.141         3.895           1.799
##      age_cat20-21 age_cat22-23 age_cat24-25 age_cat26-29 age_cat30-34
## [1,]        0.919        0.946        0.880        0.799        0.687
## [2,]        1.088        1.174        1.082        1.117        0.998
##      age_cat35-49 age_cat50-64 age_cat65+ sexualityBisexual sexualityLes/Gay
## [1,]        0.552        0.356      0.305             1.761            1.246
## [2,]        0.829        0.511      0.494             1.500            1.046
##      maleMale marstdivorced marstseparated marstwidowed lst_alc_use212>1month
## [1,]    0.818         1.281          1.649        1.559                 0.951
## [2,]    0.886         0.893          1.921        1.745                 0.815
##      lst_alc_use2last 30days dep_year21
## [1,]                   0.851     10.521
## [2,]                   0.782      9.765

Multinomial Regresion Work

mfit<-vglm(hopelesslst30days~race_eth+educ+age_cat+sexuality+male+marst+lst_alc_use2+dep_year2,
           family=multinomial(refLevel = 1),
           data=NSDUH_2019,
           weights=analwtc/mean(analwtc, na.rm=T))

summary(mfit)
## 
## Call:
## vglm(formula = hopelesslst30days ~ race_eth + educ + age_cat + 
##     sexuality + male + marst + lst_alc_use2 + dep_year2, family = multinomial(refLevel = 1), 
##     data = NSDUH_2019, weights = analwtc/mean(analwtc, na.rm = T))
## 
## Coefficients: 
##                           Estimate Std. Error z value Pr(>|z|)    
## (Intercept):1             -0.57222    0.09621  -5.948 2.72e-09 ***
## (Intercept):2             -3.75524    0.19195 -19.564  < 2e-16 ***
## race_ethasian:1           -0.13296    0.21263  -0.625 0.531759    
## race_ethasian:2            0.33098    0.42351   0.782 0.434494    
## race_ethblack:1           -0.23768    0.04057  -5.858 4.69e-09 ***
## race_ethblack:2           -0.02054    0.08647  -0.237 0.812274    
## race_ethhispanic:1         0.25532    0.08387   3.044 0.002332 ** 
## race_ethhispanic:2         0.26149    0.16715   1.564 0.117708    
## race_ethmult_race:1        0.15305    0.05804   2.637 0.008367 ** 
## race_ethmult_race:2        0.40133    0.14389   2.789 0.005284 ** 
## race_ethother:1           -0.20697    0.17134  -1.208 0.227069    
## race_ethother:2            0.38064    0.29195   1.304 0.192304    
## educassociates:1           0.13762    0.04575   3.008 0.002630 ** 
## educassociates:2           0.47986    0.12441   3.857 0.000115 ***
## educhighschool:1           0.07897    0.03650   2.163 0.030505 *  
## educhighschool:2           1.18661    0.09067  13.087  < 2e-16 ***
## educLssThnHgh:1            0.25087    0.05318   4.717 2.39e-06 ***
## educLssThnHgh:2            1.48168    0.11446  12.944  < 2e-16 ***
## educsomeCollege:1          0.09360    0.03598   2.601 0.009287 ** 
## educsomeCollege:2          0.64132    0.09361   6.851 7.33e-12 ***
## age_cat20-21:1            -0.09807    0.10478  -0.936 0.349313    
## age_cat20-21:2             0.01823    0.17586   0.104 0.917457    
## age_cat22-23:1            -0.06912    0.10153  -0.681 0.496043    
## age_cat22-23:2             0.10486    0.17466   0.600 0.548266    
## age_cat24-25:1            -0.13368    0.10111  -1.322 0.186122    
## age_cat24-25:2            -0.01831    0.17887  -0.102 0.918475    
## age_cat26-29:1            -0.23491    0.08998  -2.611 0.009040 ** 
## age_cat26-29:2            -0.06831    0.15776  -0.433 0.665004    
## age_cat30-34:1            -0.38152    0.09004  -4.237 2.26e-05 ***
## age_cat30-34:2            -0.25764    0.16180  -1.592 0.111301    
## age_cat35-49:1            -0.59433    0.08653  -6.868 6.50e-12 ***
## age_cat35-49:2            -0.54692    0.15488  -3.531 0.000414 ***
## age_cat50-64:1            -1.01250    0.08789 -11.520  < 2e-16 ***
## age_cat50-64:2            -1.19667    0.16424  -7.286 3.19e-13 ***
## age_cat65+:1              -1.17908    0.09098 -12.960  < 2e-16 ***
## age_cat65+:2              -1.22459    0.17862  -6.856 7.09e-12 ***
## sexualityBisexual:1        0.52747    0.06581   8.015 1.11e-15 ***
## sexualityBisexual:2        0.80830    0.10592   7.631 2.33e-14 ***
## sexualityLes/Gay:1         0.22658    0.08847   2.561 0.010436 *  
## sexualityLes/Gay:2         0.17702    0.18697   0.947 0.343742    
## maleMale:1                -0.19934    0.02642  -7.545 4.52e-14 ***
## maleMale:2                -0.21504    0.06249  -3.441 0.000579 ***
## marstdivorced:1            0.27791    0.06133   4.532 5.85e-06 ***
## marstdivorced:2           -0.05243    0.18297  -0.287 0.774466    
## marstseparated:1           0.46551    0.03735  12.465  < 2e-16 ***
## marstseparated:2           0.85347    0.08645   9.872  < 2e-16 ***
## marstwidowed:1             0.41962    0.03895  10.774  < 2e-16 ***
## marstwidowed:2             0.71943    0.09309   7.728 1.09e-14 ***
## lst_alc_use212>1month:1   -0.03468    0.04413  -0.786 0.431944    
## lst_alc_use212>1month:2   -0.20578    0.09843  -2.091 0.036558 *  
## lst_alc_use2last 30days:1 -0.14737    0.03535  -4.169 3.05e-05 ***
## lst_alc_use2last 30days:2 -0.30220    0.08087  -3.737 0.000186 ***
## dep_year21:1               2.16479    0.05089  42.542  < 2e-16 ***
## dep_year21:2               3.50315    0.07293  48.032  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Names of linear predictors: log(mu[,2]/mu[,1]), log(mu[,3]/mu[,1])
## 
## Residual deviance: 45390.6 on 59280 degrees of freedom
## 
## Log-likelihood: -22695.3 on 59280 degrees of freedom
## 
## Number of Fisher scoring iterations: 7 
## 
## Warning: Hauck-Donner effect detected in the following estimate(s):
## '(Intercept):2'
## 
## 
## Reference group is level  1  of the response
round(exp(coef(mfit)), 3)
##             (Intercept):1             (Intercept):2           race_ethasian:1 
##                     0.564                     0.023                     0.876 
##           race_ethasian:2           race_ethblack:1           race_ethblack:2 
##                     1.392                     0.788                     0.980 
##        race_ethhispanic:1        race_ethhispanic:2       race_ethmult_race:1 
##                     1.291                     1.299                     1.165 
##       race_ethmult_race:2           race_ethother:1           race_ethother:2 
##                     1.494                     0.813                     1.463 
##          educassociates:1          educassociates:2          educhighschool:1 
##                     1.148                     1.616                     1.082 
##          educhighschool:2           educLssThnHgh:1           educLssThnHgh:2 
##                     3.276                     1.285                     4.400 
##         educsomeCollege:1         educsomeCollege:2            age_cat20-21:1 
##                     1.098                     1.899                     0.907 
##            age_cat20-21:2            age_cat22-23:1            age_cat22-23:2 
##                     1.018                     0.933                     1.111 
##            age_cat24-25:1            age_cat24-25:2            age_cat26-29:1 
##                     0.875                     0.982                     0.791 
##            age_cat26-29:2            age_cat30-34:1            age_cat30-34:2 
##                     0.934                     0.683                     0.773 
##            age_cat35-49:1            age_cat35-49:2            age_cat50-64:1 
##                     0.552                     0.579                     0.363 
##            age_cat50-64:2              age_cat65+:1              age_cat65+:2 
##                     0.302                     0.308                     0.294 
##       sexualityBisexual:1       sexualityBisexual:2        sexualityLes/Gay:1 
##                     1.695                     2.244                     1.254 
##        sexualityLes/Gay:2                maleMale:1                maleMale:2 
##                     1.194                     0.819                     0.807 
##           marstdivorced:1           marstdivorced:2          marstseparated:1 
##                     1.320                     0.949                     1.593 
##          marstseparated:2            marstwidowed:1            marstwidowed:2 
##                     2.348                     1.521                     2.053 
##   lst_alc_use212>1month:1   lst_alc_use212>1month:2 lst_alc_use2last 30days:1 
##                     0.966                     0.814                     0.863 
## lst_alc_use2last 30days:2              dep_year21:1              dep_year21:2 
##                     0.739                     8.713                    33.220
round(exp(confint(mfit)), 3)
##                            2.5 % 97.5 %
## (Intercept):1              0.467  0.681
## (Intercept):2              0.016  0.034
## race_ethasian:1            0.577  1.328
## race_ethasian:2            0.607  3.193
## race_ethblack:1            0.728  0.854
## race_ethblack:2            0.827  1.161
## race_ethhispanic:1         1.095  1.521
## race_ethhispanic:2         0.936  1.802
## race_ethmult_race:1        1.040  1.306
## race_ethmult_race:2        1.127  1.981
## race_ethother:1            0.581  1.138
## race_ethother:2            0.826  2.593
## educassociates:1           1.049  1.255
## educassociates:2           1.266  2.062
## educhighschool:1           1.007  1.162
## educhighschool:2           2.743  3.913
## educLssThnHgh:1            1.158  1.426
## educLssThnHgh:2            3.516  5.507
## educsomeCollege:1          1.023  1.178
## educsomeCollege:2          1.581  2.281
## age_cat20-21:1             0.738  1.113
## age_cat20-21:2             0.721  1.438
## age_cat22-23:1             0.765  1.139
## age_cat22-23:2             0.789  1.564
## age_cat24-25:1             0.718  1.067
## age_cat24-25:2             0.692  1.394
## age_cat26-29:1             0.663  0.943
## age_cat26-29:2             0.686  1.272
## age_cat30-34:1             0.572  0.815
## age_cat30-34:2             0.563  1.061
## age_cat35-49:1             0.466  0.654
## age_cat35-49:2             0.427  0.784
## age_cat50-64:1             0.306  0.432
## age_cat50-64:2             0.219  0.417
## age_cat65+:1               0.257  0.368
## age_cat65+:2               0.207  0.417
## sexualityBisexual:1        1.490  1.928
## sexualityBisexual:2        1.823  2.762
## sexualityLes/Gay:1         1.055  1.492
## sexualityLes/Gay:2         0.827  1.722
## maleMale:1                 0.778  0.863
## maleMale:2                 0.714  0.912
## marstdivorced:1            1.171  1.489
## marstdivorced:2            0.663  1.358
## marstseparated:1           1.480  1.714
## marstseparated:2           1.982  2.781
## marstwidowed:1             1.410  1.642
## marstwidowed:2             1.711  2.464
## lst_alc_use212>1month:1    0.886  1.053
## lst_alc_use212>1month:2    0.671  0.987
## lst_alc_use2last 30days:1  0.805  0.925
## lst_alc_use2last 30days:2  0.631  0.866
## dep_year21:1               7.886  9.627
## dep_year21:2              28.795 38.325

Fitted Probabilities

dat<-expand.grid(race_eth=levels(NSDUH_2019$race_eth),
                 educ=levels(NSDUH_2019$educ),
                 male=levels(NSDUH_2019$male),
                 age_cat=levels(NSDUH_2019$age_cat),
                marst=levels(NSDUH_2019$marst),
                sexuality=levels(NSDUH_2019$sexuality),
                lst_alc_use2=levels(NSDUH_2019$lst_alc_use2),
                dep_year2=levels(NSDUH_2019$dep_year2))
fitm<-predict(mfit, newdat=dat,type="response")
dat<-cbind(dat, round(fitm,3))
head(dat, n=20)
##     race_eth       educ   male age_cat   marst    sexuality lst_alc_use2
## 1      white    colgrad Female   18-19 married Heterosexual    >12months
## 2      asian    colgrad Female   18-19 married Heterosexual    >12months
## 3      black    colgrad Female   18-19 married Heterosexual    >12months
## 4   hispanic    colgrad Female   18-19 married Heterosexual    >12months
## 5  mult_race    colgrad Female   18-19 married Heterosexual    >12months
## 6      other    colgrad Female   18-19 married Heterosexual    >12months
## 7      white associates Female   18-19 married Heterosexual    >12months
## 8      asian associates Female   18-19 married Heterosexual    >12months
## 9      black associates Female   18-19 married Heterosexual    >12months
## 10  hispanic associates Female   18-19 married Heterosexual    >12months
## 11 mult_race associates Female   18-19 married Heterosexual    >12months
## 12     other associates Female   18-19 married Heterosexual    >12months
## 13     white highschool Female   18-19 married Heterosexual    >12months
## 14     asian highschool Female   18-19 married Heterosexual    >12months
## 15     black highschool Female   18-19 married Heterosexual    >12months
## 16  hispanic highschool Female   18-19 married Heterosexual    >12months
## 17 mult_race highschool Female   18-19 married Heterosexual    >12months
## 18     other highschool Female   18-19 married Heterosexual    >12months
## 19     white  LssThnHgh Female   18-19 married Heterosexual    >12months
## 20     asian  LssThnHgh Female   18-19 married Heterosexual    >12months
##    dep_year2     1     2     3
## 1          0 0.630 0.355 0.015
## 2          0 0.655 0.324 0.021
## 3          0 0.681 0.303 0.016
## 4          0 0.569 0.414 0.017
## 5          0 0.591 0.389 0.021
## 6          0 0.670 0.307 0.023
## 7          0 0.593 0.384 0.022
## 8          0 0.617 0.350 0.032
## 9          0 0.646 0.330 0.024
## 10         0 0.531 0.443 0.026
## 11         0 0.552 0.417 0.031
## 12         0 0.632 0.333 0.035
## 13         0 0.593 0.362 0.045
## 14         0 0.609 0.326 0.065
## 15         0 0.642 0.309 0.048
## 16         0 0.530 0.418 0.053
## 17         0 0.548 0.390 0.063
## 18         0 0.622 0.309 0.070
## 19         0 0.547 0.397 0.056
## 20         0 0.562 0.357 0.081

Non-Proportional and proportional Model work

## Proportional
fit.vgam<-vglm(as.ordered(hopelesslst30days)~race_eth+educ+age_cat+sexuality+male+marst+lst_alc_use2+dep_year2,
               NSDUH_2019, weights =analwtc/mean(analwtc, na.rm=T),
               family=cumulative(parallel = T, reverse = T))  #<-parallel = T == proportional odds
summary(fit.vgam)
## 
## Call:
## vglm(formula = as.ordered(hopelesslst30days) ~ race_eth + educ + 
##     age_cat + sexuality + male + marst + lst_alc_use2 + dep_year2, 
##     family = cumulative(parallel = T, reverse = T), data = NSDUH_2019, 
##     weights = analwtc/mean(analwtc, na.rm = T))
## 
## Coefficients: 
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept):1           -0.59532    0.08758  -6.798 1.06e-11 ***
## (Intercept):2           -3.46474    0.09226 -37.556  < 2e-16 ***
## race_ethasian           -0.03944    0.19682  -0.200 0.841163    
## race_ethblack           -0.18957    0.03786  -5.007 5.52e-07 ***
## race_ethhispanic         0.20189    0.07735   2.610 0.009051 ** 
## race_ethmult_race        0.17677    0.05523   3.201 0.001372 ** 
## race_ethother           -0.04832    0.15335  -0.315 0.752687    
## educassociates           0.15844    0.04391   3.608 0.000309 ***
## educhighschool           0.24401    0.03433   7.108 1.18e-12 ***
## educLssThnHgh            0.43087    0.04921   8.756  < 2e-16 ***
## educsomeCollege          0.15207    0.03423   4.442 8.92e-06 ***
## age_cat20-21            -0.06813    0.09303  -0.732 0.463985    
## age_cat22-23            -0.05236    0.09062  -0.578 0.563419    
## age_cat24-25            -0.12523    0.09075  -1.380 0.167576    
## age_cat26-29            -0.19005    0.08070  -2.355 0.018516 *  
## age_cat30-34            -0.33137    0.08100  -4.091 4.30e-05 ***
## age_cat35-49            -0.54241    0.07777  -6.975 3.06e-12 ***
## age_cat50-64            -0.99122    0.07932 -12.497  < 2e-16 ***
## age_cat65+              -1.14614    0.08252 -13.889  < 2e-16 ***
## sexualityBisexual        0.48741    0.05676   8.587  < 2e-16 ***
## sexualityLes/Gay         0.20272    0.08216   2.467 0.013610 *  
## maleMale                -0.19185    0.02503  -7.665 1.79e-14 ***
## marstdivorced            0.23308    0.05922   3.936 8.30e-05 ***
## marstseparated           0.50798    0.03513  14.460  < 2e-16 ***
## marstwidowed             0.44614    0.03680  12.122  < 2e-16 ***
## lst_alc_use212>1month   -0.06224    0.04148  -1.500 0.133503    
## lst_alc_use2last 30days -0.16664    0.03333  -4.999 5.77e-07 ***
## dep_year21               2.28639    0.04129  55.375  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Names of linear predictors: logitlink(P[Y>=2]), logitlink(P[Y>=3])
## 
## Residual deviance: 45651.43 on 59306 degrees of freedom
## 
## Log-likelihood: -22825.71 on 59306 degrees of freedom
## 
## Number of Fisher scoring iterations: 5 
## 
## No Hauck-Donner effect found in any of the estimates
## 
## 
## Exponentiated coefficients:
##           race_ethasian           race_ethblack        race_ethhispanic 
##               0.9613245               0.8273160               1.2237142 
##       race_ethmult_race           race_ethother          educassociates 
##               1.1933600               0.9528283               1.1716772 
##          educhighschool           educLssThnHgh         educsomeCollege 
##               1.2763630               1.5385958               1.1642359 
##            age_cat20-21            age_cat22-23            age_cat24-25 
##               0.9341433               0.9489874               0.8822924 
##            age_cat26-29            age_cat30-34            age_cat35-49 
##               0.8269151               0.7179383               0.5813463 
##            age_cat50-64              age_cat65+       sexualityBisexual 
##               0.3711250               0.3178623               1.6280974 
##        sexualityLes/Gay                maleMale           marstdivorced 
##               1.2247296               0.8254338               1.2624886 
##          marstseparated            marstwidowed   lst_alc_use212>1month 
##               1.6619338               1.5622666               0.9396572 
## lst_alc_use2last 30days              dep_year21 
##               0.8465074               9.8393820
#Non-proportional odds
fit.vgam2<-vglm(as.ordered(hopelesslst30days)~race_eth+educ+age_cat+sexuality+male+marst+lst_alc_use2+dep_year2,
                NSDUH_2019,
                weights =analwtc/mean(analwtc, na.rm=T),
                family=cumulative(parallel = F, reverse = T))  #<-parallel = F == Nonproportional odds
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1

## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1

## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in Deviance.categorical.data.vgam(mu = mu, y = y, w = w, residuals =
## residuals, : fitted values close to 0 or 1
## Warning in slot(family, "validparams")(eta, y, extra = extra): It seems that
## the nonparallelism assumption has resulted in intersecting linear/additive
## predictors. Try propodds() or fitting a partial nonproportional odds model or
## choosing some other link function, etc.
## Warning in vglm.fitter(x = x, y = y, w = w, offset = offset, Xm2 = Xm2, :
## iterations terminated because half-step sizes are very small
## Warning in vglm.fitter(x = x, y = y, w = w, offset = offset, Xm2 = Xm2, : some
## quantities such as z, residuals, SEs may be inaccurate due to convergence at a
## half-step
## Warning in log(prob): NaNs produced
summary(fit.vgam2)
## Warning in matrix.power(wz, M = M, power = 0.5, fast = TRUE): Some weight
## matrices have negative eigenvalues. They will be assigned NAs
## 
## Call:
## vglm(formula = as.ordered(hopelesslst30days) ~ race_eth + educ + 
##     age_cat + sexuality + male + marst + lst_alc_use2 + dep_year2, 
##     family = cumulative(parallel = F, reverse = T), data = NSDUH_2019, 
##     weights = analwtc/mean(analwtc, na.rm = T))
## 
## Coefficients: 
##                             Estimate Std. Error    z value Pr(>|z|)    
## (Intercept):1             -4.476e-01  1.282e-02 -3.492e+01   <2e-16 ***
## (Intercept):2             -3.313e+00  3.121e-02 -1.062e+02   <2e-16 ***
## race_ethasian:1           -6.378e-02  2.112e-07 -3.020e+05   <2e-16 ***
## race_ethasian:2            3.953e-01  2.646e-07  1.494e+06   <2e-16 ***
## race_ethblack:1           -1.832e-01  5.130e-08 -3.572e+06   <2e-16 ***
## race_ethblack:2            6.748e-02  6.239e-08  1.082e+06   <2e-16 ***
## race_ethhispanic:1         2.324e-01  1.270e-07  1.830e+06   <2e-16 ***
## race_ethhispanic:2         1.523e-01  1.645e-07  9.259e+05   <2e-16 ***
## race_ethmult_race:1        1.416e-01  1.453e-07  9.745e+05   <2e-16 ***
## race_ethmult_race:2        2.386e-01  1.861e-07  1.282e+06   <2e-16 ***
## race_ethother:1           -1.053e-01  1.051e-07 -1.001e+06   <2e-16 ***
## race_ethother:2            7.172e-01  1.447e-07  4.955e+06   <2e-16 ***
## educassociates:1           1.254e-01  2.733e-07  4.590e+05   <2e-16 ***
## educassociates:2           1.977e-01  3.441e-07  5.745e+05   <2e-16 ***
## educhighschool:1           1.389e-01  2.170e-07  6.402e+05   <2e-16 ***
## educhighschool:2           9.107e-01  2.566e-07  3.549e+06   <2e-16 ***
## educLssThnHgh:1            3.068e-01  2.267e-07  1.353e+06   <2e-16 ***
## educLssThnHgh:2            1.209e+00  2.614e-07  4.625e+06   <2e-16 ***
## educsomeCollege:1          1.012e-01  2.096e-07  4.826e+05   <2e-16 ***
## educsomeCollege:2          3.166e-01  2.607e-07  1.215e+06   <2e-16 ***
## age_cat20-21:1            -1.041e-01  8.068e-08 -1.291e+06   <2e-16 ***
## age_cat20-21:2            -2.792e-02  1.161e-07 -2.404e+05   <2e-16 ***
## age_cat22-23:1            -8.407e-02  1.083e-07 -7.760e+05   <2e-16 ***
## age_cat22-23:2            -1.750e-01  1.624e-07 -1.078e+06   <2e-16 ***
## age_cat24-25:1            -1.663e-01  1.231e-07 -1.351e+06   <2e-16 ***
## age_cat24-25:2            -4.010e-01  1.754e-07 -2.286e+06   <2e-16 ***
## age_cat26-29:1            -2.628e-01  9.850e-08 -2.668e+06   <2e-16 ***
## age_cat26-29:2            -4.465e-01  1.383e-07 -3.227e+06   <2e-16 ***
## age_cat30-34:1            -4.320e-01  1.049e-07 -4.119e+06   <2e-16 ***
## age_cat30-34:2            -6.437e-01  1.415e-07 -4.549e+06   <2e-16 ***
## age_cat35-49:1            -6.350e-01  8.995e-08 -7.059e+06   <2e-16 ***
## age_cat35-49:2            -7.957e-01  1.238e-07 -6.429e+06   <2e-16 ***
## age_cat50-64:1            -9.879e-01  1.000e-07 -9.875e+06   <2e-16 ***
## age_cat50-64:2            -1.098e+00  1.281e-07 -8.570e+06   <2e-16 ***
## age_cat65+:1              -1.097e+00  9.803e-08 -1.119e+07   <2e-16 ***
## age_cat65+:2              -1.018e+00  1.325e-07 -7.685e+06   <2e-16 ***
## sexualityBisexual:1        5.010e-01  8.745e-08  5.729e+06   <2e-16 ***
## sexualityBisexual:2        1.215e+00  8.463e-08  1.435e+07   <2e-16 ***
## sexualityLes/Gay:1         1.950e-01  1.545e-07  1.263e+06   <2e-16 ***
## sexualityLes/Gay:2         5.372e-03  2.021e-07  2.657e+04   <2e-16 ***
## maleMale:1                -1.636e-01  4.183e-08 -3.911e+06   <2e-16 ***
## maleMale:2                -9.250e-02  5.339e-08 -1.733e+06   <2e-16 ***
## marstdivorced:1            1.789e-01  1.850e-07  9.671e+05   <2e-16 ***
## marstdivorced:2           -2.395e-01  2.121e-07 -1.129e+06   <2e-16 ***
## marstseparated:1           4.455e-01  6.730e-08  6.619e+06   <2e-16 ***
## marstseparated:2           5.311e-01  7.419e-08  7.159e+06   <2e-16 ***
## marstwidowed:1             3.663e-01  7.261e-08  5.045e+06   <2e-16 ***
## marstwidowed:2             3.435e-01  9.074e-08  3.786e+06   <2e-16 ***
## lst_alc_use212>1month:1   -3.777e-02  6.718e-08 -5.622e+05   <2e-16 ***
## lst_alc_use212>1month:2   -1.267e-01  8.513e-08 -1.488e+06   <2e-16 ***
## lst_alc_use2last 30days:1 -1.322e-01  4.874e-08 -2.713e+06   <2e-16 ***
## lst_alc_use2last 30days:2 -1.623e-01  6.174e-08 -2.629e+06   <2e-16 ***
## dep_year21:1               2.329e+00  1.282e-02  1.817e+02   <2e-16 ***
## dep_year21:2               4.905e+00  3.121e-02  1.572e+02   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Names of linear predictors: logitlink(P[Y>=2]), logitlink(P[Y>=3])
## 
## Residual deviance: 594353.6 on 59280 degrees of freedom
## 
## Log-likelihood: NA on 59280 degrees of freedom
## 
## Number of Fisher scoring iterations: 2 
## 
## Warning: Hauck-Donner effect detected in the following estimate(s):
## '(Intercept):2', 'race_ethhispanic:1', 'race_ethmult_race:1', 'educassociates:1', 'educhighschool:1', 'educLssThnHgh:1', 'educsomeCollege:1', 'age_cat20-21:2', 'age_cat22-23:2', 'age_cat24-25:2', 'age_cat26-29:2', 'age_cat30-34:2', 'age_cat35-49:2', 'age_cat50-64:2', 'age_cat65+:2', 'sexualityBisexual:1', 'sexualityLes/Gay:1', 'maleMale:2', 'marstdivorced:1', 'marstdivorced:2', 'marstseparated:1', 'marstwidowed:1', 'lst_alc_use212>1month:2', 'lst_alc_use2last 30days:2', 'dep_year21:1'
## 
## 
## Exponentiated coefficients:
##           race_ethasian:1           race_ethasian:2           race_ethblack:1 
##                 0.9382134                 1.4848387                 0.8325749 
##           race_ethblack:2        race_ethhispanic:1        race_ethhispanic:2 
##                 1.0698126                 1.2616653                 1.1645589 
##       race_ethmult_race:1       race_ethmult_race:2           race_ethother:1 
##                 1.1521373                 1.2694318                 0.9000594 
##           race_ethother:2          educassociates:1          educassociates:2 
##                 2.0487340                 1.1336577                 1.2185531 
##          educhighschool:1          educhighschool:2           educLssThnHgh:1 
##                 1.1490023                 2.4860187                 1.3590598 
##           educLssThnHgh:2         educsomeCollege:1         educsomeCollege:2 
##                 3.3501447                 1.1064579                 1.3725138 
##            age_cat20-21:1            age_cat20-21:2            age_cat22-23:1 
##                 0.9011048                 0.9724697                 0.9193700 
##            age_cat22-23:2            age_cat24-25:1            age_cat24-25:2 
##                 0.8394280                 0.8468018                 0.6696217 
##            age_cat26-29:1            age_cat26-29:2            age_cat30-34:1 
##                 0.7689212                 0.6398798                 0.6492363 
##            age_cat30-34:2            age_cat35-49:1            age_cat35-49:2 
##                 0.5253237                 0.5299591                 0.4512819 
##            age_cat50-64:1            age_cat50-64:2              age_cat65+:1 
##                 0.3723636                 0.3335572                 0.3340076 
##              age_cat65+:2       sexualityBisexual:1       sexualityBisexual:2 
##                 0.3613094                 1.6504152                 3.3699845 
##        sexualityLes/Gay:1        sexualityLes/Gay:2                maleMale:1 
##                 1.2153596                 1.0053864                 0.8490645 
##                maleMale:2           marstdivorced:1           marstdivorced:2 
##                 0.9116501                 1.1959512                 0.7870223 
##          marstseparated:1          marstseparated:2            marstwidowed:1 
##                 1.5612556                 1.7007923                 1.4423905 
##            marstwidowed:2   lst_alc_use212>1month:1   lst_alc_use212>1month:2 
##                 1.4099432                 0.9629342                 0.8810107 
## lst_alc_use2last 30days:1 lst_alc_use2last 30days:2              dep_year21:1 
##                 0.8761585                 0.8501458                10.2701789 
##              dep_year21:2 
##               135.0039272
AIC(fit.vgam)
## [1] 45707.43
AIC(fit.vgam2)
## [1] NaN
attach(dat)
fitted.ord<-round(predict(fit.vgam, newdat=dat[,1:3], type="response"), 3)

dat<-cbind(dat, fitted.ord)

names(dat)<-c(names(dat)[1:3], "mp1", "mp2", "mp3", "op1", "op2", "op3")

head(dat, n=20)
##     race_eth       educ   male   mp1     mp2          mp3       op1 op2   op3
## 1      white    colgrad Female 18-19 married Heterosexual >12months   0 0.630
## 2      asian    colgrad Female 18-19 married Heterosexual >12months   0 0.655
## 3      black    colgrad Female 18-19 married Heterosexual >12months   0 0.681
## 4   hispanic    colgrad Female 18-19 married Heterosexual >12months   0 0.569
## 5  mult_race    colgrad Female 18-19 married Heterosexual >12months   0 0.591
## 6      other    colgrad Female 18-19 married Heterosexual >12months   0 0.670
## 7      white associates Female 18-19 married Heterosexual >12months   0 0.593
## 8      asian associates Female 18-19 married Heterosexual >12months   0 0.617
## 9      black associates Female 18-19 married Heterosexual >12months   0 0.646
## 10  hispanic associates Female 18-19 married Heterosexual >12months   0 0.531
## 11 mult_race associates Female 18-19 married Heterosexual >12months   0 0.552
## 12     other associates Female 18-19 married Heterosexual >12months   0 0.632
## 13     white highschool Female 18-19 married Heterosexual >12months   0 0.593
## 14     asian highschool Female 18-19 married Heterosexual >12months   0 0.609
## 15     black highschool Female 18-19 married Heterosexual >12months   0 0.642
## 16  hispanic highschool Female 18-19 married Heterosexual >12months   0 0.530
## 17 mult_race highschool Female 18-19 married Heterosexual >12months   0 0.548
## 18     other highschool Female 18-19 married Heterosexual >12months   0 0.622
## 19     white  LssThnHgh Female 18-19 married Heterosexual >12months   0 0.547
## 20     asian  LssThnHgh Female 18-19 married Heterosexual >12months   0 0.562
##       NA    NA    NA    NA    NA
## 1  0.355 0.015 0.645 0.325 0.030
## 2  0.324 0.021 0.654 0.317 0.029
## 3  0.303 0.016 0.687 0.288 0.025
## 4  0.414 0.017 0.597 0.366 0.037
## 5  0.389 0.021 0.603 0.361 0.036
## 6  0.307 0.023 0.656 0.315 0.029
## 7  0.384 0.022 0.608 0.357 0.035
## 8  0.350 0.032 0.617 0.349 0.034
## 9  0.330 0.024 0.652 0.319 0.029
## 10 0.443 0.026 0.558 0.399 0.043
## 11 0.417 0.031 0.565 0.393 0.042
## 12 0.333 0.035 0.619 0.347 0.034
## 13 0.362 0.045 0.587 0.375 0.038
## 14 0.326 0.065 0.596 0.367 0.037
## 15 0.309 0.048 0.632 0.336 0.032
## 16 0.418 0.053 0.537 0.416 0.047
## 17 0.390 0.063 0.544 0.411 0.045
## 18 0.309 0.070 0.599 0.365 0.037
## 19 0.397 0.056 0.541 0.413 0.046
## 20 0.357 0.081 0.551 0.405 0.044

##Are the assumptions of the proportional odds model met? How did you evaluate the proportional odds assumption? ## A sensitivity analysis was performed on both the ordinal and multinomial regression models. The ordinal model did not have much consistent proportionality, as much variation in the observed odds ratios occured in the education variables, and many of the racial based variables. The most consistent of these appeared to be the mental health, alcohol, and sexuality based variables. Comparing this to the multinomial model the multinomial one has less variation in each of the categorical variables used ti eplain the outcome variable. It appears to be more consistent but less of a fit for the model.

##Which model (Proportional odds or Multinomial) fits the data better? how did you decide upon this? ## Comparing the AIC reveals the multinomial model to be the best fit for my data. Comparing both the proportional odds model and multinomial indicates an aic of 45707.43, vs 45498.6, the lower of the two indicates a better fit.

AIC(fit.vgam)
## [1] 45707.43
AIC(mfit)
## [1] 45498.6

##For the best fitting model from part c,Describe the results of your model, and present output from the model in terms of odds ratios and confidence intervals for all model parameters in a table. ##The results of the multinomial model indicate several interesting items. Compared to whites, blacks were less likely to report not being hopeless, vs sometimes feeling hopeless. Hispanics were more likely than whites to feel hopeless some of the time vs not at all, with multi racial people feeling hopeless some of the time, compared to not often, and again feeling hopeless most of the time vs some of the time compared to whites. Compared to people with College degrees, those without them including all other other educational attainment levels (less than high school, high school, some college, associates) reported feeling hopeles some of the time vs none of the time, and more often feeling hopeless than sometimes feeling hopeless. Education means lower levels of feeling hopeless not often. People aged 26-29 and 30-34 were less likely to feel hopeless some of the time vs not often at all compared to people aged 18-19. Compared to people aged 18-19 , people aged 35-49, 50-64, and 65+ were are less likely to report feeling hopeless some of the time vs none of the time, and most of the time compared to some of the time. Homosexuals compared to heterosexuals are more likely to indicate feeling hopeless some of the time as opposed to none of time. However bisexuals compared to heterosexuals are both more likely to report feeling hopeless some of the time compared to not often feeling hopeless, and feeling hopeless most fo the time vs some of the time. Males are less likely to report feeling hopeless some of the time vs none of the time, and most of the time vs some of the time than their female counterparts. Compared to married individuals, divorced people are more likely to report feeling hopeless some of the time vs none of the time. WIth these results being similar for separated and widowed individuals, but these people also are more likely to report feeling hopeless most of the time vs some of the time. Alcohol usage within the last year is less likely to result in feeling hopless most of the time vs some of the time, however if alcohol was used within the last thrity days the individual will be less likely to report feeling hopless some of the time vs not often, and most of the time vs some of the time. Finally depressive episodes within the last year is more likely to result in feelings of hopelessness some of the time, vs none of the time, and most of the time compared to some.

mfit<-vglm(hopelesslst30days~race_eth+educ+age_cat+sexuality+male+marst+lst_alc_use2+dep_year2,
           family=multinomial(refLevel = 1),
           data=NSDUH_2019,
           weights=analwtc/mean(analwtc, na.rm=T))

summary(mfit)
## 
## Call:
## vglm(formula = hopelesslst30days ~ race_eth + educ + age_cat + 
##     sexuality + male + marst + lst_alc_use2 + dep_year2, family = multinomial(refLevel = 1), 
##     data = NSDUH_2019, weights = analwtc/mean(analwtc, na.rm = T))
## 
## Coefficients: 
##                           Estimate Std. Error z value Pr(>|z|)    
## (Intercept):1             -0.57222    0.09621  -5.948 2.72e-09 ***
## (Intercept):2             -3.75524    0.19195 -19.564  < 2e-16 ***
## race_ethasian:1           -0.13296    0.21263  -0.625 0.531759    
## race_ethasian:2            0.33098    0.42351   0.782 0.434494    
## race_ethblack:1           -0.23768    0.04057  -5.858 4.69e-09 ***
## race_ethblack:2           -0.02054    0.08647  -0.237 0.812274    
## race_ethhispanic:1         0.25532    0.08387   3.044 0.002332 ** 
## race_ethhispanic:2         0.26149    0.16715   1.564 0.117708    
## race_ethmult_race:1        0.15305    0.05804   2.637 0.008367 ** 
## race_ethmult_race:2        0.40133    0.14389   2.789 0.005284 ** 
## race_ethother:1           -0.20697    0.17134  -1.208 0.227069    
## race_ethother:2            0.38064    0.29195   1.304 0.192304    
## educassociates:1           0.13762    0.04575   3.008 0.002630 ** 
## educassociates:2           0.47986    0.12441   3.857 0.000115 ***
## educhighschool:1           0.07897    0.03650   2.163 0.030505 *  
## educhighschool:2           1.18661    0.09067  13.087  < 2e-16 ***
## educLssThnHgh:1            0.25087    0.05318   4.717 2.39e-06 ***
## educLssThnHgh:2            1.48168    0.11446  12.944  < 2e-16 ***
## educsomeCollege:1          0.09360    0.03598   2.601 0.009287 ** 
## educsomeCollege:2          0.64132    0.09361   6.851 7.33e-12 ***
## age_cat20-21:1            -0.09807    0.10478  -0.936 0.349313    
## age_cat20-21:2             0.01823    0.17586   0.104 0.917457    
## age_cat22-23:1            -0.06912    0.10153  -0.681 0.496043    
## age_cat22-23:2             0.10486    0.17466   0.600 0.548266    
## age_cat24-25:1            -0.13368    0.10111  -1.322 0.186122    
## age_cat24-25:2            -0.01831    0.17887  -0.102 0.918475    
## age_cat26-29:1            -0.23491    0.08998  -2.611 0.009040 ** 
## age_cat26-29:2            -0.06831    0.15776  -0.433 0.665004    
## age_cat30-34:1            -0.38152    0.09004  -4.237 2.26e-05 ***
## age_cat30-34:2            -0.25764    0.16180  -1.592 0.111301    
## age_cat35-49:1            -0.59433    0.08653  -6.868 6.50e-12 ***
## age_cat35-49:2            -0.54692    0.15488  -3.531 0.000414 ***
## age_cat50-64:1            -1.01250    0.08789 -11.520  < 2e-16 ***
## age_cat50-64:2            -1.19667    0.16424  -7.286 3.19e-13 ***
## age_cat65+:1              -1.17908    0.09098 -12.960  < 2e-16 ***
## age_cat65+:2              -1.22459    0.17862  -6.856 7.09e-12 ***
## sexualityBisexual:1        0.52747    0.06581   8.015 1.11e-15 ***
## sexualityBisexual:2        0.80830    0.10592   7.631 2.33e-14 ***
## sexualityLes/Gay:1         0.22658    0.08847   2.561 0.010436 *  
## sexualityLes/Gay:2         0.17702    0.18697   0.947 0.343742    
## maleMale:1                -0.19934    0.02642  -7.545 4.52e-14 ***
## maleMale:2                -0.21504    0.06249  -3.441 0.000579 ***
## marstdivorced:1            0.27791    0.06133   4.532 5.85e-06 ***
## marstdivorced:2           -0.05243    0.18297  -0.287 0.774466    
## marstseparated:1           0.46551    0.03735  12.465  < 2e-16 ***
## marstseparated:2           0.85347    0.08645   9.872  < 2e-16 ***
## marstwidowed:1             0.41962    0.03895  10.774  < 2e-16 ***
## marstwidowed:2             0.71943    0.09309   7.728 1.09e-14 ***
## lst_alc_use212>1month:1   -0.03468    0.04413  -0.786 0.431944    
## lst_alc_use212>1month:2   -0.20578    0.09843  -2.091 0.036558 *  
## lst_alc_use2last 30days:1 -0.14737    0.03535  -4.169 3.05e-05 ***
## lst_alc_use2last 30days:2 -0.30220    0.08087  -3.737 0.000186 ***
## dep_year21:1               2.16479    0.05089  42.542  < 2e-16 ***
## dep_year21:2               3.50315    0.07293  48.032  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Names of linear predictors: log(mu[,2]/mu[,1]), log(mu[,3]/mu[,1])
## 
## Residual deviance: 45390.6 on 59280 degrees of freedom
## 
## Log-likelihood: -22695.3 on 59280 degrees of freedom
## 
## Number of Fisher scoring iterations: 7 
## 
## Warning: Hauck-Donner effect detected in the following estimate(s):
## '(Intercept):2'
## 
## 
## Reference group is level  1  of the response
round(exp(coef(mfit)), 3)
##             (Intercept):1             (Intercept):2           race_ethasian:1 
##                     0.564                     0.023                     0.876 
##           race_ethasian:2           race_ethblack:1           race_ethblack:2 
##                     1.392                     0.788                     0.980 
##        race_ethhispanic:1        race_ethhispanic:2       race_ethmult_race:1 
##                     1.291                     1.299                     1.165 
##       race_ethmult_race:2           race_ethother:1           race_ethother:2 
##                     1.494                     0.813                     1.463 
##          educassociates:1          educassociates:2          educhighschool:1 
##                     1.148                     1.616                     1.082 
##          educhighschool:2           educLssThnHgh:1           educLssThnHgh:2 
##                     3.276                     1.285                     4.400 
##         educsomeCollege:1         educsomeCollege:2            age_cat20-21:1 
##                     1.098                     1.899                     0.907 
##            age_cat20-21:2            age_cat22-23:1            age_cat22-23:2 
##                     1.018                     0.933                     1.111 
##            age_cat24-25:1            age_cat24-25:2            age_cat26-29:1 
##                     0.875                     0.982                     0.791 
##            age_cat26-29:2            age_cat30-34:1            age_cat30-34:2 
##                     0.934                     0.683                     0.773 
##            age_cat35-49:1            age_cat35-49:2            age_cat50-64:1 
##                     0.552                     0.579                     0.363 
##            age_cat50-64:2              age_cat65+:1              age_cat65+:2 
##                     0.302                     0.308                     0.294 
##       sexualityBisexual:1       sexualityBisexual:2        sexualityLes/Gay:1 
##                     1.695                     2.244                     1.254 
##        sexualityLes/Gay:2                maleMale:1                maleMale:2 
##                     1.194                     0.819                     0.807 
##           marstdivorced:1           marstdivorced:2          marstseparated:1 
##                     1.320                     0.949                     1.593 
##          marstseparated:2            marstwidowed:1            marstwidowed:2 
##                     2.348                     1.521                     2.053 
##   lst_alc_use212>1month:1   lst_alc_use212>1month:2 lst_alc_use2last 30days:1 
##                     0.966                     0.814                     0.863 
## lst_alc_use2last 30days:2              dep_year21:1              dep_year21:2 
##                     0.739                     8.713                    33.220
round(exp(confint(mfit)), 3)
##                            2.5 % 97.5 %
## (Intercept):1              0.467  0.681
## (Intercept):2              0.016  0.034
## race_ethasian:1            0.577  1.328
## race_ethasian:2            0.607  3.193
## race_ethblack:1            0.728  0.854
## race_ethblack:2            0.827  1.161
## race_ethhispanic:1         1.095  1.521
## race_ethhispanic:2         0.936  1.802
## race_ethmult_race:1        1.040  1.306
## race_ethmult_race:2        1.127  1.981
## race_ethother:1            0.581  1.138
## race_ethother:2            0.826  2.593
## educassociates:1           1.049  1.255
## educassociates:2           1.266  2.062
## educhighschool:1           1.007  1.162
## educhighschool:2           2.743  3.913
## educLssThnHgh:1            1.158  1.426
## educLssThnHgh:2            3.516  5.507
## educsomeCollege:1          1.023  1.178
## educsomeCollege:2          1.581  2.281
## age_cat20-21:1             0.738  1.113
## age_cat20-21:2             0.721  1.438
## age_cat22-23:1             0.765  1.139
## age_cat22-23:2             0.789  1.564
## age_cat24-25:1             0.718  1.067
## age_cat24-25:2             0.692  1.394
## age_cat26-29:1             0.663  0.943
## age_cat26-29:2             0.686  1.272
## age_cat30-34:1             0.572  0.815
## age_cat30-34:2             0.563  1.061
## age_cat35-49:1             0.466  0.654
## age_cat35-49:2             0.427  0.784
## age_cat50-64:1             0.306  0.432
## age_cat50-64:2             0.219  0.417
## age_cat65+:1               0.257  0.368
## age_cat65+:2               0.207  0.417
## sexualityBisexual:1        1.490  1.928
## sexualityBisexual:2        1.823  2.762
## sexualityLes/Gay:1         1.055  1.492
## sexualityLes/Gay:2         0.827  1.722
## maleMale:1                 0.778  0.863
## maleMale:2                 0.714  0.912
## marstdivorced:1            1.171  1.489
## marstdivorced:2            0.663  1.358
## marstseparated:1           1.480  1.714
## marstseparated:2           1.982  2.781
## marstwidowed:1             1.410  1.642
## marstwidowed:2             1.711  2.464
## lst_alc_use212>1month:1    0.886  1.053
## lst_alc_use212>1month:2    0.671  0.987
## lst_alc_use2last 30days:1  0.805  0.925
## lst_alc_use2last 30days:2  0.631  0.866
## dep_year21:1               7.886  9.627
## dep_year21:2              28.795 38.325