|
Dear R/S users,
This is related to my previous mail.
I came with the bellow way to obtain the SE’s
Although, it give the correct estimates I think it is very inefficient,
especially when we have 4 year of data and aprox. 70 variables to analyze per year.
If someone came up with something better, please let me know.
Sincerely,
Jose
library(foreign)
library(survey)
BRFSS.06<-read.spss("C:/Datos2006AsmaC1.sav",
use.value.labels = TRUE, to.data.frame = TRUE)
BRFSS.06.svy<-svydesign(ids=~1, strata=~STSTR, weights=~FINALWT,
data="">
BRFSS.06.svy<-as.svydesign2(BRFSS.06.svy)
B.06.Sset.svy<- subset (BRFSS.06.svy, ASTHMA2==1)
Nw.06<-sum(BRFSS.06$FINALWT)
ANow.mean<-((svytotal(~factor(ASTHNOW), design=B.06.Sset,
na.rm=TRUE, estimate.
#estimate.only do not work
ANow.mean<-matrix(ANow.mean)[1,1]
ANow.SE<-1.96 *((SE(svytotal(~factor(ASTHNOW),
design=BRFSS.06.svy, na.rm=TRUE))/Nw.06)*100)
ANow.SE<-matrix(ANow.SE)[1,1]
ANow.UCL<- ANow.mean + ANow.SE
ANow.LCL<- ANow.mean - ANow.SE
ANow.Est<- c(ANow.UCL, ANow.mean, ANow.LCL)
Est.names<-c("upper Conf Lim", "Prevalence", "Lower Conf Lim")
names(ANow.Est)<-Est.names
ANow.Est
####### For tabulated analysis ###
AsthNow.Sset<-data.frame(subset(BRFSS.06, ASTHNOW==1,
select=c(ASTHMA2, ASTHNOW, SEX, EDUCAG, INCOMG,
STSTR, FINALWT)))
dim(AsthNow.Sset)
attach(AsthNow.Sset)
AsthNow.Sset.svy<-svydesign(ids=~1, strata=~STSTR,
weights=~FINALWT, data="">
AsthNow.Sset.svy<-as.svydesign2(AsthNow.Sset.svy)
####Mean####
sum(ASTHNOW * FINALWT)
svytotal(~ASTHNOW, design=AsthNow.Sset.svy)/Nw.06
AstH.Sex.Num<-svyby(~ASTHNOW, ~factor(SEX), AsthNow.Sset.svy,
svytotal)
AstH.Sex.Num<-data.frame(AstH.Sex.Num)
AstH.Sex.Num[1:2,2]
ANow.Sex.Mean<-(AstH.Sex.Num[1:2,2]/tapply(BRFSS.06$FINALWT,
BRFSS.06$SEX, sum))*100
####SE###
AstH.Sex.Num[1:2,3]
AstH.Sex.LCL<-ANow.Sex.Mean - 1.96 *(AstH.Sex.Num[1:2,3]/tapply(BRFSS.06$FINALWT, BRFSS.06$SEX, sum)* 100)
AstH.Sex.UCL<-ANow.Sex.Mean + 1.96 *(AstH.Sex.Num[1:2,3]/tapply(BRFSS.06$FINALWT, BRFSS.06$SEX, sum)* 100)
AsthNow.Male.Est<- c(AstH.Sex.LCL[1], ANow.Sex.Mean[1], AstH.Sex.UCL[1])
names(AsthNow.Male.Est)<-Est.names
AsthNow.Female.Est<- c(AstH.Sex.LCL[2], ANow.Sex.Mean[2], AstH.Sex.UCL[2])
names(AsthNow.Female.Est)<-Est.names
AsthNow.Male.Est
AsthNow.Female.Est
Use video conversation to talk face-to-face with Windows Live Messenger. Get started.
|