Analysis of health and drug use
Red text is commands and blue is results
> summary(ex5.des)
Independent Sampling design
svydesign(id = ~1, weights = ~WEIGHT, data = ex5)
Probabilities:
Min. 1st Qu. Median Mean 3rd Qu. Max.
0.3282 0.7363 0.8910 0.8962 1.0960 1.3170
Data variables:
[1] "AGE" "WEIGHT" "EMPL" "QUAL" "GENHELF" "Q85A" "Q85B" "SACC" "SINC" "GWEIGHT" "LIVING"
> #
> # now we can use the svymean commands to get first the weighted proportions
> # for categorical variables
> svymean(~GENHELF,ex5.des,deff=T)
mean SE DEff
GENHELFexcellent 0.1866598 0.0207625 1.0250
GENHELFvery good 0.3776246 0.0266726 1.0928
GENHELFgood 0.3435661 0.0263967 1.1153
GENHELFfair 0.0700424 0.0143276 1.1377
GENHELFpoor 0.0221072 0.0094428 1.4890
> svymean(~Q85A,ex5.des,deff=T)
mean SE DEff
Q85Anever used NA NA NA
Q85Atried once or twice NA NA NA
Q85Ause daily NA NA NA
Q85Ause weekly NA NA NA
Q85Aused in last month NA NA NA
Q85Aused more than a month ago NA NA NA
> #
> # this fails for missing values so need to get subset of design without missing values
> #
> svymean(~Q85A,subset(ex5.des,!is.na(Q85A)),deff=T)
mean SE DEff
Q85Anever used 0.5831140 0.0274052 1.1061
Q85Atried once or twice 0.2497045 0.0240431 1.1046
Q85Ause daily 0.0340624 0.0101163 1.1135
Q85Ause weekly 0.0240867 0.0088028 1.1801
Q85Aused in last month 0.0277697 0.0084588 0.9488
Q85Aused more than a month ago 0.0812627 0.0159420 1.2187
> #
> # and expressing this as percentages
> #
> round(print(svymean(~Q85A,subset(ex5.des,!is.na(Q85A))))*100,1)
mean SE
Q85Anever used 58.3 2.7
Q85Atried once or twice 25.0 2.4
Q85Ause daily 3.4 1.0
Q85Ause weekly 2.4 0.9
Q85Aused in last month 2.8 0.8
Q85Aused more than a month ago 8.1 1.6
>
> # now the means of some continuous variables
> #
> svymean(~as.numeric(GENHELF),ex5.des,deff=T) # this treats it as a score not a category
mean SE DEff
as.numeric(GENHELF) 2.363313 0.052741 1.1445
> svymean(~SINC,ex5.des,deff=T)
mean SE DEff
SINC 18.7117 0.6946 1.3158
> svymean(~SACC,ex5.des,deff=T)
mean SE DEff
SACC -0.067572 0.032883 1.0413
> summary(svyglm(as.numeric(GENHELF)~CANSCORE,ex5.des))
Call:
svyglm(formula = as.numeric(GENHELF) ~ CANSCORE, design = ex5.des)
Survey design:
update.survey.design(ex5.des, AMPSCORE = AMPSCORE)
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.28376 0.05666 40.304 < 2e-16 ***
CANSCORE 0.45467 0.14982 3.035 0.00258 **
---
Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1
(Dispersion parameter for gaussian family taken to be 0.803459)
Number of Fisher Scoring iterations: 2
> #
> # Now make a table to understand results
> #
> can.table<-svytable(~GENHELF+Q85A,ex5.des)
> can.tprint<-round(sweep(can.table,2,apply(can.table,2,sum),"/")*100) # calculating column percents
> can.tprint<-rbind(can.tprint,table(ex5$Q85A)) # add a row below for base numbers
> #
> # at any stage type tyhe name of an object to see what it is
> #
> dimnames(can.tprint)[[1]][6]<-'BASE'
> can.tprint
never once/twice daily weekly last month >month ago
excellent 22 19 6 9 0 8
very good 40 35 41 30 17 42
good 29 37 53 42 54 48
fair 7 8 0 19 20 2
poor 2 0 0 0 8 0
BASE 210 90 12 8 11 27
>
>