# NOT RUN {
data(datenKapitel02)
schueler <- datenKapitel02$schueler
schule <- datenKapitel02$schule
set.seed(20150506)
# }
# NOT RUN {
## -------------------------------------------------------------
## Abschnitt 4.1: Stratifizierung
## -------------------------------------------------------------
# -------------------------------------------------------------
# Abschnitt 4.1, Listing 1
# Information in Strata
strata <- aggregate(schule[,"NSchueler", drop = FALSE],
by=schule[,"stratum", drop = FALSE], sum)
colnames(strata)[2] <- "NSchuelerStratum"
## -------------------------------------------------------------
## Abschnitt 4.2: Schulenziehung
## -------------------------------------------------------------
# -------------------------------------------------------------
# Abschnitt 4.2, Listing 1
# Dummyvariable Klassenziehung
schule$Klassenziehung <- 0
schule[which(schule$NKlassen>3), "Klassenziehung"] <- 1
# -------------------------------------------------------------
# Abschnitt 4.2, Listing 2
# erwarteter Beitrag zur Stichprobe pro Schule
schule$NSchueler.erw <- schule$NSchueler
ind <- which(schule$Klassenziehung == 1)
schule[ind, "NSchueler.erw"] <-
schule[ind, "NSchueler"]/schule[ind, "NKlassen"]*3
# -------------------------------------------------------------
# Abschnitt 4.2, Listing 3
# relativer Anteil Sch<U+00FC>ler pro Schule
temp <- merge(schule[, c("SKZ","stratum","NSchueler")],
strata[, c("stratum","NSchuelerStratum")])
schule$AnteilSchueler <-
temp$NSchueler/temp$NSchuelerStratum
# mittlere Anzahl von Sch<U+00FC>lern pro Schule
strata$"NSchueler/Schule.erw" <-
rowsum(apply(schule, 1, function(x)
x["NSchueler.erw"]*x["AnteilSchueler"]), schule$stratum)
# -------------------------------------------------------------
# Abschnitt 4.2, Listing 4
# Bestimmung Anzahl zu ziehender Schulen
strata$Schulen.zu.ziehen <-
round(2500/strata[,"NSchueler/Schule.erw"])
# -------------------------------------------------------------
# Abschnitt 4.2, Listing 5
# Schulenliste nach Stratum und Groesse ordnen
schule <-
schule[order(schule$stratum, schule$NSchueler),]
# -------------------------------------------------------------
# Abschnitt 4.2, Listing 6
# Berechnung Sampling-Intervall
strata$Samp.Int <-
strata$NSchuelerStratum/strata$Schulen.zu.ziehen
# -------------------------------------------------------------
# Abschnitt 4.2, Listing 7
# Startwerte bestimmen
strata$Startwert <-
sapply(ceiling(strata$Samp.Int), sample, size = 1)
# -------------------------------------------------------------
# Abschnitt 4.2, Listing 8
# Sch<U+00FC>ler-Tickets
tickets <- sapply(1:4, function(x)
trunc(0:(strata[strata$stratum==x,"Schulen.zu.ziehen"]-1)
* strata[strata$stratum==x, "Samp.Int"] +
strata$Startwert[x]))
# -------------------------------------------------------------
# Abschnitt 4.2, Listing 9
# kummulierte Sch<U+00FC>leranzahl pro Stratum berechnen
schule$NSchuelerKum <-
unlist(sapply(1:4, function(x)
cumsum(schule[schule$stratum==x, "NSchueler"])))
# -------------------------------------------------------------
# Abschnitt 4.2, Listing 10
# Schulen ziehen
schule$SInSamp <- 0
for(s in 1:4) {
NSchuelerKumStrat <-
schule[schule$stratum==s, "NSchuelerKum"]
inds <- sapply(tickets[[s]], function(x)
setdiff(which(NSchuelerKumStrat <= x),
which(NSchuelerKumStrat[-1] <= x)))
schule[schule$stratum==s, "SInSamp"][inds] <- 1 }
# -------------------------------------------------------------
# Abschnitt 4.2, Listing 11
# Berechnung Ziehungswahrscheinlichkeit Schule
temp <- merge(schule[, c("stratum", "AnteilSchueler")],
strata[, c("stratum", "Schulen.zu.ziehen")])
schule$Z.Wsk.Schule <-
temp$AnteilSchueler*temp$Schulen.zu.ziehen
## -------------------------------------------------------------
## Abschnitt 4.3: Klassenziehung
## -------------------------------------------------------------
# -------------------------------------------------------------
# Abschnitt 4.3, Listing 1
### Klassenziehung (Alternative 2)
schukla <- unique(merge(
schule[, c("SKZ","NKlassen", "Klassenziehung",
"Z.Wsk.Schule", "SInSamp")],
schueler[, c("SKZ", "idclass")], by="SKZ"))
schukla$KlInSamp <- 0
for(skz in unique(schukla[schukla$SInSamp==1,"SKZ"])) {
temp <- schukla[schukla$SKZ==skz, "idclass"]
schukla[schukla$idclass%in%temp[sample.int(
min(3, length(temp)))], "KlInSamp"] <- 1 }
# -------------------------------------------------------------
# Abschnitt 4.3, Listing 2
# Ziehungswahrscheinlichkeit Klasse
schukla$Z.Wsk.Klasse <- ((1 - schukla$Klassenziehung) * 1 +
schukla$Klassenziehung * 3 / schukla$NKlassen)
## -------------------------------------------------------------
## Abschnitt 4.4: Gewichtung
## -------------------------------------------------------------
# -------------------------------------------------------------
# Abschnitt 4.4, Listing 1
### Gewichte
schueler <- merge(schueler, schukla[, c("idclass", "KlInSamp", "Z.Wsk.Schule",
"Z.Wsk.Klasse")],
by="idclass", all.x=T)
# Ziehungswahrscheinlichkeiten Schueler
schueler$Z.Wsk.Schueler <-
schueler$Z.Wsk.Schule * schueler$Z.Wsk.Klasse
# -------------------------------------------------------------
# Abschnitt 4.4, Listing 2
schueler <- schueler[schueler$KlInSamp==1,]
# Nonresponse Adjustment
temp <- merge(schueler[, c("idclass", "Z.Wsk.Schueler")],
aggregate(schueler$teilnahme,
by=list(schueler$idclass),
function(x) sum(x)/length(x)),
by.x="idclass", by.y="Group.1")
# -------------------------------------------------------------
# Abschnitt 4.4, Listing 3
# Sch<U+00FC>lergewichte
schueler$studwgt <- 1/temp$x/temp$Z.Wsk.Schueler
# -------------------------------------------------------------
# Abschnitt 4.4, Listing 4
# Normierung
Normierung <- strata$NSchuelerStratum /
rowsum(schueler$studwgt * schueler$teilnahme,
group = schueler$Stratum)
schueler$NormStudwgt <-
schueler$studwgt * Normierung[schueler$Stratum]
## -------------------------------------------------------------
## Abschnitt 5.3: Jackknife-Repeated-Replication
## -------------------------------------------------------------
# -------------------------------------------------------------
# Abschnitt 5.3, Listing 1
### Erg<U+00E4>nzung zum Buch: Hilfsfunktion zones.within.stratum
zones.within.stratum <- function(offset,n.str) {
maxzone <- offset-1+floor(n.str/2)
zones <- sort(rep(offset:maxzone,2))
if (n.str %% 2 == 1) zones <- c(zones,maxzone)
return(zones) }
### Ende der Erg<U+00E4>nzung
# Sortieren der Schulliste (explizite und implizite Strata)
schule <- schule[schule$SInSamp==1,]
schule <- schule[order(schule$stratum,-schule$NSchueler),]
# Unterteilung in Pseudostrata
cnt.strata <- length(unique(schule$stratum))
offset <- 1
jkzones.vect <- integer()
for (i in 1:cnt.strata) {
n.str <- table(schule$stratum)[i]
jkzones.vect <-
c(jkzones.vect,zones.within.stratum(offset,n.str))
offset <- max(jkzones.vect)+1 }
schule$jkzone <- jkzones.vect
# Zuf<U+00E4>llige Auswahl von Schulen mit Gewicht 0
schule$jkrep <- 1
cnt.zones <- max(schule$jkzone)
jkrep.rows.null <- integer()
for (i in 1:cnt.zones) {
rows.zone <- which(schule$jkzone==i)
### Erg<U+00E4>nzung zum Buch: Fallunterscheidung je nach Anzahl Schulen in der Zone
if (length(rows.zone)==2) jkrep.rows.null <-
c(jkrep.rows.null,sample(rows.zone,size=1))
else {
num.null <- sample(1:2,size=1)
jkrep.rows.null <-
c(jkrep.rows.null,sample(rows.zone,size=num.null))
} }
schule[jkrep.rows.null,]$jkrep <- 0
# -------------------------------------------------------------
# Abschnitt 5.3, Listing 2
# <U+00DC>bertragung auf Sch<U+00FC>lerebene
schueler <-
merge(schueler,schule[,c("SKZ","jkzone","jkrep")],all.x=TRUE)
# Schleife zur Generierung von Replicate Weights
for (i in 1:cnt.zones) {
in.zone <- as.numeric(schueler$jkzone==i)
schueler[paste0("w_fstr",i)] <- # vgl. Formel 5
in.zone * schueler$jkrep * schueler$NormStudwgt * 2 +
(1-in.zone) * schueler$NormStudwgt }
# -------------------------------------------------------------
# Abschnitt 5.3, Listing 3
# Sch<U+00E4>tzung mittels Gesamtgewicht
n.female <- sum(schueler[schueler$female==1,]$NormStudwgt)
perc.female <- n.female / sum(schueler$NormStudwgt)
# wiederholte Berechnung und Varianz
var.jrr = 0
for (i in 1:cnt.zones) {
n.female.rep <-
sum(schueler[schueler$female==1,paste0("w_fstr",i)])
perc.female.rep <-
n.female.rep / sum(schueler[paste0("w_fstr",i)])
var.jrr <- # vgl. Formel 6
var.jrr + (perc.female.rep - perc.female) ^ 2.0 }
# }
# NOT RUN {
# }
# NOT RUN {
<!-- % end dontrun -->
# }
Run the code above in your browser using DataLab