handout7

library(SDaA)
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

par(mfrow=c(1,2))
data(htpop)
popecdf <- ecdf(htpop$height)
plot(popecdf, do.points = FALSE, ylab = "F(y)", xlab = "Height Value, y")

minht <- min(htpop$height)
breaks <- c(minht - 1, seq(from = minht, to = max(htpop$height), by = 1))
hist(htpop$height, ylab = "f(y)", breaks = breaks,
xlab = "Height Value, y", freq = FALSE)

par(mfrow=c(1,3))
minht <- min(htpop$height)
breaks <- c(minht - 1, seq(from = minht, to = max(htpop$height), by = 1))
hist(htpop$height, ylab = "f(y)", breaks = breaks,
xlab = "Height Value, y", freq = FALSE)


data(htsrs)
hist(htsrs$height, ylab = "Relative Frequency", xlab = "Height (cm)",
freq = FALSE)

data(htstrat)
hist(htstrat$height, ylab = "Relative Frequency",
xlab = "Height (cm)", freq = FALSE)

par(mfrow=c(1,2))
minht <- min(htstrat$height)
breaks <- c(minht - 1, seq(from = minht, to = max(htstrat$height),
by = 1))
hist(htstrat$height, ylab = expression(hat(f)(y)),
breaks = breaks, xlab = "Height Value, y", freq = FALSE)

stratecdf <- ecdf(htstrat$height)
plot(stratecdf, do.points = FALSE, ylab = expression(hat(F)(y)),
xlab = "Height Value, y")

data(htpop)
mean(htpop$height)

## [1] 168.616

quantile(htpop$height,probs=c(25,50,90)/100)

## 25% 50% 90% 
## 160 168 184

var(htpop$height)

## [1] 124.5038

data(htsrs)
mean(htsrs$height)

## [1] 168.94

quantile(htsrs$height,probs=c(25,50,90)/100)

## 25% 50% 90% 
## 160 169 184

var(htsrs$height)

## [1] 122.6396

data(htstrat)
mean(htstrat$height)

## [1] 164.645

quantile(htstrat$height,probs=c(25,50,90)/100)

## 25% 50% 90% 
## 157 163 178

var(htstrat$height)

## [1] 93.40601

str<-c(rep(1,160),rep(2,40))
strdesign<-svydesign(id = ~1, strata = ~str, fpc = c(rep(1000,160),rep(1000,40)), 
data = htstrat)
svymean(~height,strdesign)

##          mean    SE
## height 169.02 0.717

svyquantile(~height,strdesign,c(.25,.50,.90))

##        0.25 0.5 0.9
## height  161 168 182

svyvar(htstrat$height,strdesign)

##      variance     SE
## [1,]   117.35 15.733

data(syc)
par(mfrow=c(1,2))
hist(syc$age, main="Sample unweighted",col="purple",
freq = FALSE)

oo <- options(survey.lonely.psu = "certainty")
sycdesign <- svydesign(ids = ~psu, strata = ~stratum,
data = syc, weights = ~finalwt)
svyhist(~age, sycdesign, main="Sample weighted",col="purple")

options(oo)


sycdesign <- svydesign(ids = ~psu, strata = ~stratum,
data = syc, weights = ~finalwt)
oo <- options(survey.lonely.psu = "certainty")
svyboxplot(age ~ factor(stratum), design = sycdesign)
options(oo)


par(mfrow=c(1,1))

library(ggplot2)

## 
## Attaching package: 'ggplot2'

## The following object is masked from 'package:SDaA':
## 
##     seals

p <- ggplot(syc, aes(x = factor(stratum), y = factor(age)))
g <- p + stat_sum(aes(group = 1, weight = finalwt,
size = ..n..))
print(g)

oo <- options(survey.lonely.psu = "certainty")
sycstrat5 <- subset(sycdesign, stratum == 5)
svyboxplot(age ~ factor(psu), design = sycstrat5)

options(oo)

sycstrat5df <- subset(syc, stratum == 5)
p <- ggplot(sycstrat5df, aes(x = factor(psu), y = factor(age)))
g <- p + stat_sum(aes(group = 1, weight = finalwt,
size = ..n..))
print(g)