library(car)
## Loading required package: carData
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:graphics':
##
## dotchart
library(sjPlot)
## Warning: package 'sjPlot' was built under R version 4.0.4
## Registered S3 methods overwritten by 'lme4':
## method from
## cooks.distance.influence.merMod car
## influence.merMod car
## dfbeta.influence.merMod car
## dfbetas.influence.merMod car
library(ggplot2)
load("~/ipumstor3/data.RData")
summary(data)
## STRATA PSU PERWEIGHT SAMPWEIGHT
## Min. :6001 Min. :1.000 Min. : 0 Min. : 0
## 1st Qu.:6081 1st Qu.:1.000 1st Qu.: 2093 1st Qu.: 0
## Median :6158 Median :2.000 Median : 3356 Median : 0
## Mean :6156 Mean :1.507 Mean : 3380 Mean : 3380
## 3rd Qu.:6238 3rd Qu.:2.000 3rd Qu.: 4394 3rd Qu.: 5466
## Max. :6300 Max. :2.000 Max. :28756 Max. :104188
## FWEIGHT SEX MARSTAT ALCAMT
## Min. : 897 Min. :1.000 Min. : 0.0 Min. : 0.000
## 1st Qu.: 1917 1st Qu.:1.000 1st Qu.:10.0 1st Qu.: 0.000
## Median : 3123 Median :2.000 Median :10.0 Median : 0.000
## Mean : 3184 Mean :1.516 Mean :20.8 Mean : 5.678
## 3rd Qu.: 4184 3rd Qu.:2.000 3rd Qu.:50.0 3rd Qu.: 0.000
## Max. :26845 Max. :2.000 Max. :99.0 Max. :999.000
## CIGSDAY
## Min. : 1.00
## 1st Qu.:96.00
## Median :96.00
## Mean :91.27
## 3rd Qu.:96.00
## Max. :99.00
#Using IPUMs data we use the count variable cigarettes per day since we also used alcohol consumed per day we did not need a offset variable
#research question, does marital status relate to smoking cigarettes
#results show that the best model is the poisson model fit2 because it has the lowest AIC
data$gender<-recode(data$SEX,
recodes="1=1; 2=0; else=NA")
#1=married, 2=widowed, 3=divorced, 4=seperated, 5=NM
data$mar<-recode(data$MARSTAT,
recodes="10:13 =1; 20=2; 30=2; 40=4; 50=5; else=NA")
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
sub<-data%>%
select(gender, mar, MARSTAT,
SEX,CIGSDAY, ALCAMT, FWEIGHT, PERWEIGHT, PSU,
SAMPWEIGHT, STRATA) %>%
filter( complete.cases(.))
#First we tell R our survey design
options(survey.lonely.psu = "adjust")
des<-svydesign(ids=~1, strata=~STRATA, weights=~SAMPWEIGHT, data =sub )
#OR THE BRFSS, R GAVE ME A WARNING AND I NEEDED TO ADD:
#YOU MAY NOT NEED TO DO THIS!!!!
#First we tell R our survey design
options(survey.lonely.psu = "adjust")
des<-svydesign(ids=~1, strata=~STRATA,
weights=~SAMPWEIGHT,
data = data[is.na(data$SAMPWEIGHT)==F,])
svyhist(~CIGSDAY, des)

svyby(~CIGSDAY, ~mar+gender, des, svymean, na.rm=T)
## mar gender CIGSDAY se
## 1.0 1 0 84.78112 0.4192275
## 2.0 2 0 79.48747 0.6204091
## 4.0 4 0 74.34465 1.7711725
## 5.0 5 0 81.80322 0.5862962
## 1.1 1 1 83.22563 0.4470258
## 2.1 2 1 71.93251 0.9536398
## 4.1 4 1 64.18817 2.5132948
## 5.1 5 1 78.05711 0.6367085
svyby(~CIGSDAY, ~ALCAMT, des, svymean, na.rm=T)
## ALCAMT CIGSDAY se
## 0 0 90.55462 0.1769934
## 10 10 85.47080 0.4300409
## 20 20 79.78102 0.5350991
## 30 30 70.97386 0.9690191
## 40 40 67.08518 1.3840793
## 50 50 61.28640 2.0790308
## 60 60 56.09178 1.9782904
## 70 70 59.32912 4.7994932
## 80 80 55.55627 4.1955653
## 90 90 51.58078 7.7145131
## 100 100 57.75151 4.3975446
## 110 110 62.69831 28.8588833
## 120 120 58.67694 4.3756821
## 130 130 96.00000 0.0000000
## 140 140 69.92419 19.7054872
## 150 150 37.56220 9.6902262
## 160 160 96.00000 0.0000000
## 180 180 25.32328 12.7048950
## 200 200 52.81000 17.5101385
## 240 240 84.97549 7.5211834
## 250 250 40.00000 0.0000000
## 300 300 82.81882 11.4038980
## 350 350 20.00000 0.0000000
## 600 600 1.00000 0.0000000
## 870 870 10.00000 0.0000000
## 950 950 65.49671 27.0041024
## 997 997 96.00000 0.0000000
## 998 998 64.06579 28.3919526
## 999 999 72.37294 4.7130137
fit1<-svyglm(CIGSDAY~factor(mar)+factor(gender), design=des, family=poisson)
summary(fit1)
##
## Call:
## svyglm(formula = CIGSDAY ~ factor(mar) + factor(gender), design = des,
## family = poisson)
##
## Survey design:
## svydesign(ids = ~1, strata = ~STRATA, weights = ~SAMPWEIGHT,
## data = data[is.na(data$SAMPWEIGHT) == F, ])
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.451252 0.004442 1002.098 < 2e-16 ***
## factor(mar)2 -0.095191 0.007793 -12.215 < 2e-16 ***
## factor(mar)4 -0.182648 0.021326 -8.565 < 2e-16 ***
## factor(mar)5 -0.050006 0.006567 -7.615 2.71e-14 ***
## factor(gender)1 -0.040979 0.005618 -7.294 3.08e-13 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 12.88573)
##
## Number of Fisher Scoring iterations: 4
round(exp(summary(fit1)$coef[-1,1]), 3)
## factor(mar)2 factor(mar)4 factor(mar)5 factor(gender)1
## 0.909 0.833 0.951 0.960
fit2<-glm(CIGSDAY~factor(mar)+factor(gender), data=data, family=poisson)
summary(fit2)
##
## Call:
## glm(formula = CIGSDAY ~ factor(mar) + factor(gender), family = poisson,
## data = data)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -13.2115 0.3301 0.4611 0.7552 1.5290
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.5304653 0.0006775 6687.40 <2e-16 ***
## factor(mar)2 -0.0795048 0.0011940 -66.59 <2e-16 ***
## factor(mar)4 -0.1063431 0.0026095 -40.75 <2e-16 ***
## factor(mar)5 -0.0306387 0.0008925 -34.33 <2e-16 ***
## factor(gender)1 -0.0135545 0.0007986 -16.97 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 599581 on 71017 degrees of freedom
## Residual deviance: 593674 on 71013 degrees of freedom
## (18958 observations deleted due to missingness)
## AIC: 1036483
##
## Number of Fisher Scoring iterations: 4
fit3<-glm(CIGSDAY~factor(mar)+factor(gender)+factor(ALCAMT), data=data, family=quasipoisson)
summary(fit3)
##
## Call:
## glm(formula = CIGSDAY ~ factor(mar) + factor(gender) + factor(ALCAMT),
## family = quasipoisson, data = data)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -13.4191 0.0413 0.2313 0.5941 9.4359
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.558095 0.001584 2877.938 < 2e-16 ***
## factor(mar)2 -0.055006 0.002704 -20.341 < 2e-16 ***
## factor(mar)4 -0.072830 0.005890 -12.364 < 2e-16 ***
## factor(mar)5 -0.019484 0.002019 -9.651 < 2e-16 ***
## factor(gender)1 0.002034 0.001811 1.123 0.261315
## factor(ALCAMT)10 -0.095045 0.003355 -28.328 < 2e-16 ***
## factor(ALCAMT)20 -0.175898 0.003689 -47.685 < 2e-16 ***
## factor(ALCAMT)30 -0.271161 0.005579 -48.600 < 2e-16 ***
## factor(ALCAMT)40 -0.353632 0.008154 -43.367 < 2e-16 ***
## factor(ALCAMT)50 -0.397382 0.011363 -34.972 < 2e-16 ***
## factor(ALCAMT)60 -0.473496 0.011744 -40.319 < 2e-16 ***
## factor(ALCAMT)70 -0.447303 0.028606 -15.637 < 2e-16 ***
## factor(ALCAMT)80 -0.518704 0.025462 -20.372 < 2e-16 ***
## factor(ALCAMT)90 -0.511140 0.049415 -10.344 < 2e-16 ***
## factor(ALCAMT)100 -0.507966 0.025511 -19.912 < 2e-16 ***
## factor(ALCAMT)110 -0.568319 0.218891 -2.596 0.009424 **
## factor(ALCAMT)120 -0.437940 0.026016 -16.833 < 2e-16 ***
## factor(ALCAMT)130 0.004219 0.230007 0.018 0.985365
## factor(ALCAMT)140 -0.340621 0.193246 -1.763 0.077967 .
## factor(ALCAMT)150 -1.068695 0.076996 -13.880 < 2e-16 ***
## factor(ALCAMT)160 0.023704 0.230009 0.103 0.917920
## factor(ALCAMT)180 -1.603230 0.231211 -6.934 4.12e-12 ***
## factor(ALCAMT)200 -0.524177 0.091553 -5.725 1.04e-08 ***
## factor(ALCAMT)240 -0.292540 0.085559 -3.419 0.000628 ***
## factor(ALCAMT)250 -0.798419 0.356365 -2.240 0.025065 *
## factor(ALCAMT)300 -0.389409 0.141973 -2.743 0.006093 **
## factor(ALCAMT)350 -1.544912 0.503911 -3.066 0.002171 **
## factor(ALCAMT)600 -4.540645 2.253543 -2.015 0.043920 *
## factor(ALCAMT)870 -2.236026 0.712635 -3.138 0.001704 **
## factor(ALCAMT)950 -0.367679 0.160974 -2.284 0.022369 *
## factor(ALCAMT)997 0.041824 0.115012 0.364 0.716123
## factor(ALCAMT)998 -0.589837 0.218889 -2.695 0.007047 **
## factor(ALCAMT)999 -0.283721 0.028811 -9.848 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for quasipoisson family taken to be 5.078453)
##
## Null deviance: 599581 on 71017 degrees of freedom
## Residual deviance: 534882 on 70985 degrees of freedom
## (18958 observations deleted due to missingness)
## AIC: NA
##
## Number of Fisher Scoring iterations: 5
library(lmtest)
## Warning: package 'lmtest' was built under R version 4.0.4
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(sandwich)
## Warning: package 'sandwich' was built under R version 4.0.4
coeftest(fit2, vcov=vcovHC(fit2, type="HC1",cluster="STRATA"))
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.5304653 0.0013317 3401.9212 < 2.2e-16 ***
## factor(mar)2 -0.0795048 0.0033507 -23.7277 < 2.2e-16 ***
## factor(mar)4 -0.1063431 0.0085903 -12.3794 < 2.2e-16 ***
## factor(mar)5 -0.0306387 0.0019825 -15.4548 < 2.2e-16 ***
## factor(gender)1 -0.0135545 0.0018260 -7.4231 1.144e-13 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
library(MASS)
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
fit.nb1<-glm.nb(CIGSDAY~factor(mar),
data=sub,
weights=SAMPWEIGHT/mean(SAMPWEIGHT, na.rm=T))
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in glm.nb(CIGSDAY ~ factor(mar), data = sub, weights = SAMPWEIGHT/
## mean(SAMPWEIGHT, : alternation limit reached
fit.nb2<-glm.nb(CIGSDAY~factor(mar)+factor(gender)+factor(ALCAMT),
data=sub,
weights=SAMPWEIGHT/mean(SAMPWEIGHT, na.rm=T))
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in sqrt(1/i): NaNs produced
## Warning in glm.nb(CIGSDAY ~ factor(mar) + factor(gender) + factor(ALCAMT), :
## alternation limit reached
#clx2(fit.nb2,cluster =sub$ststr)
tests1<-coeftest(fit.nb1, vcov=vcovHC(fit.nb2, type="HC1",cluster="STRATA"))
tests<-coeftest(fit.nb2, vcov=vcovHC(fit.nb2, type="HC1",cluster="STRATA"))
library(stargazer)
stargazer(fit.nb1, fit.nb2,style="demography", type = "text", t.auto=F,p.auto=F,coef=list(tests1[, 1],tests[,1]), se =list(tests1[, 2], tests[, 2]), p=list(tests1[,4],tests[, 4]) )
##
## ---------------------------------------------
## CIGSDAY
## Model 1 Model 2
## ---------------------------------------------
## factor(mar)2 -0.089*** -0.094***
## (0.003) (0.003)
## factor(mar)4 -0.178*** -0.161***
## (0.009) (0.009)
## factor(mar)5 -0.051*** -0.037***
## (0.003) (0.003)
## factor(gender)1 0.001
## (0.002)
## factor(ALCAMT)10 -0.029***
## (0.003)
## factor(ALCAMT)20 -0.098***
## (0.003)
## factor(ALCAMT)30 -0.209***
## (0.006)
## factor(ALCAMT)40 -0.265***
## (0.009)
## factor(ALCAMT)50 -0.355***
## (0.014)
## factor(ALCAMT)60 -0.445***
## (0.015)
## factor(ALCAMT)70 -0.386***
## (0.034)
## factor(ALCAMT)80 -0.447***
## (0.032)
## factor(ALCAMT)90 -0.528***
## (0.063)
## factor(ALCAMT)100 -0.411***
## (0.032)
## factor(ALCAMT)110 -0.327
## (0.194)
## factor(ALCAMT)120 -0.395***
## (0.031)
## factor(ALCAMT)130 0.061***
## (0.002)
## factor(ALCAMT)140 -0.256*
## (0.119)
## factor(ALCAMT)150 -0.826***
## (0.107)
## factor(ALCAMT)160 0.098***
## (0.002)
## factor(ALCAMT)180 -1.253***
## (0.208)
## factor(ALCAMT)200 -0.504***
## (0.142)
## factor(ALCAMT)240 -0.011
## (0.035)
## factor(ALCAMT)250 -0.654***
## (0.009)
## factor(ALCAMT)300 -0.054
## (0.057)
## factor(ALCAMT)350 -1.471***
## (0.002)
## factor(ALCAMT)600 -4.467***
## (0.002)
## factor(ALCAMT)870 -2.163***
## (0.002)
## factor(ALCAMT)950 -0.306
## (0.169)
## factor(ALCAMT)997 0.130***
## (0.007)
## factor(ALCAMT)998 -0.343
## (0.187)
## factor(ALCAMT)999 -0.186***
## (0.027)
## Constant 4.431*** 4.502***
## (0.002) (0.002)
## N 71,018 71,018
## Log Likelihood -888,686.100 -853,615.000
## theta 878,673.400 1,503,626.000
## AIC 1,777,380.000 1,707,296.000
## ---------------------------------------------
## *p < .05; **p < .01; ***p < .001
AIC(fit.nb1)
## [1] 1777380
AIC(fit.nb2)
## [1] 1707296
AIC(fit1)
## eff.p AIC deltabar
## 1.898724e+02 1.379723e+06 4.746811e+01
AIC(fit2)
## [1] 1036483
AIC(fit3)
## [1] NA