library(survey)
data(api)
clus2_design <- svydesign(id=~dnum+snum, fpc=~fpc1+fpc2, data=apiclus2)
pop.types <- data.frame(stype=c("E","H","M"),Freq = c(4421,755,1018))
ps_design <- postStratify(clus2_design, strata=~stype,population=pop.types)
svytotal(~enroll, clus2_design,na.rm=TRUE)
svytotal(~enroll, ps_design,na.rm=TRUE)
strat_design <- svydesign(id=~1, strata=~stype, fpc=~fpc, data=apistrat)
boot_design <- as.svrepdesign(strat_design, type="bootstrap", replicates=100)
jk_design <- as.svrepdesign(strat_design, type="JKn")
brr_design <- as.svrepdesign(strat_design, type="BRR")
svytotal(~enroll, strat_design)
svytotal(~enroll, boot_design)
svytotal(~enroll, jk_design)
svytotal(~enroll, brr_design)

cdf.est<-svycdf(~enroll,strat_design)
cdf.pop<-ecdf (apipop$enroll)
cdf.samp<-ecdf (apistrat$enroll)
par(mar=c(4.1,4.1,1.1,2.1))
plot(cdf.pop,do.points=FALSE, xlab="Studenti Iscritti",ylab="Distribuzione Cumulata",main="",lwd=1)
lines(cdf.samp, do.points=FALSE,lwd=2,col="red")
lines(cdf.est[[1]] ,lwd=2,col="blue",do.points=FALSE)
legend("bottomright", lwd=c(1,2,2) ,bty="n",col=c("black","red","blue"), legend=c("Population","Weighted estimate", "Unweighted estimate") )
opar<-par(mfrow=c(1,3))
svyhist(~enroll, strat_design, main="Pesati",xlab="Studenti Iscritti",col="gray60",ylim=c(0,1.3e-3))
hist(apistrat$enroll,  main="Non Pesati",xlab="Studenti Iscritti",col="gray60",prob=TRUE,ylim=c(0,1.3e-3))
hist(apipop$enroll,  main="Popolazione",xlab="Studenti Iscritti",col="gray60",prob=TRUE,ylim=c(0,1.3e-3))
par(mfrow=c(1,2))
svyplot(api00~api99, design=strat_design, style="bubble")
svyplot(api00~api99, design=strat_design, style="transparent",pch=19)
coplot(apiclus2$api00~apiclus2$api99|apiclus2$sch.wide*apiclus2$comp.imp,pch=19)
svycoplot(api00~api99|sch.wide*comp.imp, design=clus2_design, style="transparent")