## 6.3 Replication Example

The R code below replicates multiple statistics and standard errors from the presentation by OPHI Research Officer Christop Jindra:

library(survey)
library(convey)
library(webuse)

# load the same microdata set used by Jindra in his presentation
webuse("nlsw88")

# coerce that tbl_df to a standard R data.frame
nlsw88 <- data.frame(nlsw88)

# create a collgrad column
nlsw88$collgrad <- factor( as.numeric(nlsw88$collgrad) ,
ordered = TRUE
)

# initiate a linearization-based survey design object
des_nlsw88 <- svydesign(ids = ~ 1 , data = nlsw88)
## Warning in svydesign.default(ids = ~1, data = nlsw88): No weights or
## probabilities supplied, assuming equal probability
# immediately run the convey_prep function on the survey design
des_nlsw88 <- convey_prep(des_nlsw88)

# PDF page 9
result <-
svyafc(
~ wage + collgrad + hours,
design = des_nlsw88,
cutoffs = list(4, 'college grad' , 26),
k = 1 / 3,
g = 0,
na.rm = T
)

result
##      alkire-foster     SE
## [1,]       0.36991 0.0053
attr(result, 'extra')
##        coef          SE
## H 0.8082070 0.008316807
## A 0.4576895 0.004573443
# PDF page 10
for (ks in seq(0.1, 1, .1)) {
result <-
svyafc(
~ wage + collgrad + hours,
design = des_nlsw88,
cutoffs = list(4, 'college grad' , 26),
k = ks ,
g = 0,
na.rm = T
)

print(result)

print(attr(result , 'extra'))

}
##      alkire-foster     SE
## [1,]       0.36991 0.0053
##        coef          SE
## H 0.8082070 0.008316807
## A 0.4576895 0.004573443
##      alkire-foster     SE
## [1,]       0.36991 0.0053
##        coef          SE
## H 0.8082070 0.008316807
## A 0.4576895 0.004573443
##      alkire-foster     SE
## [1,]       0.36991 0.0053
##        coef          SE
## H 0.8082070 0.008316807
## A 0.4576895 0.004573443
##      alkire-foster     SE
## [1,]       0.18659 0.0068
##        coef          SE
## H 0.2582516 0.009245464
## A 0.7225101 0.005174483
##      alkire-foster     SE
## [1,]       0.18659 0.0068
##        coef          SE
## H 0.2582516 0.009245464
## A 0.7225101 0.005174483
##      alkire-foster     SE
## [1,]       0.18659 0.0068
##        coef          SE
## H 0.2582516 0.009245464
## A 0.7225101 0.005174483
##      alkire-foster     SE
## [1,]      0.043265 0.0043
##         coef          SE
## H 0.04326494 0.004297766
## A 1.00000000 0.000000000
##      alkire-foster     SE
## [1,]      0.043265 0.0043
##         coef          SE
## H 0.04326494 0.004297766
## A 1.00000000 0.000000000
##      alkire-foster     SE
## [1,]      0.043265 0.0043
##         coef          SE
## H 0.04326494 0.004297766
## A 1.00000000 0.000000000
##      alkire-foster     SE
## [1,]      0.043265 0.0043
##         coef          SE
## H 0.04326494 0.004297766
## A 1.00000000 0.000000000
# PDF page 13
for (ks in c(0.5 , 0.75 , 1)) {
result <-
svyafc(
~ wage + collgrad + hours ,
design = des_nlsw88 ,
cutoffs = list(4, 'college grad' , 26) ,
k = ks ,
g = 0 ,
dimw = c(0.5 , 0.25 , 0.25) ,
na.rm = TRUE
)

print(result)

print(attr(result, 'extra'))

}
##      alkire-foster     SE
## [1,]       0.19135 0.0069
##        coef          SE
## H 0.2689563 0.009366804
## A 0.7114428 0.006847428
##      alkire-foster     SE
## [1,]       0.14897 0.0067
##        coef          SE
## H 0.1842105 0.008188892
## A 0.8087167 0.005216042
##      alkire-foster     SE
## [1,]      0.043265 0.0043
##         coef          SE
## H 0.04326494 0.004297766
## A 1.00000000 0.000000000
# PDF page 17
result <-
svyafc(
~ wage + collgrad + hours,
design = des_nlsw88,
cutoffs = list(4, 'college grad' , 26),
k = 1 / 3,
g = 0,
na.rm = T
)

result
##      alkire-foster     SE
## [1,]       0.36991 0.0053
attr(result, "extra")
##        coef          SE
## H 0.8082070 0.008316807
## A 0.4576895 0.004573443
svyafcdec(
~ wage + collgrad + hours,
by = ~ 1,
design = des_nlsw88,
cutoffs = list(4, 'college grad' , 26),
k = 1 / 3,
g = 0,
na.rm = T
)[2]
## $raw headcount ratio ## raw headcount SE ## wage 0.19492 0.0084 ## collgrad 0.76316 0.0090 ## hours 0.15165 0.0076 # PDF page 19 # Censored headcount ratios for (ks in seq(1 / 3, 1, 1 / 3)) { print ( svyafcdec( ~ wage + collgrad + hours, by = ~ 1, design = des_nlsw88, cutoffs = list(4, 'college grad' , 26), k = ks, g = 0, na.rm = T )[3] ) } ##$censored headcount ratio
## wage             0.19492 0.0084
## hours            0.15165 0.0076
##
## $censored headcount ratio ## cens. headcount SE ## wage 0.18421 0.0082 ## collgrad 0.24799 0.0091 ## hours 0.12756 0.0070 ## ##$censored headcount ratio
## wage            0.043265 0.0043
## hours           0.043265 0.0043
# Percentage contribution by dimensions
for (ks in seq(1 / 3, 1, 1 / 3)) {
print (
svyafcdec(
~ wage + collgrad + hours,
by = ~ 1,
design = des_nlsw88,
cutoffs = list(4, 'college grad' , 26),
k = ks,
g = 0,
na.rm = T
)[4]
)

}
## $percentual contribution per dimension ## dim. % contribution SE ## wage 0.17564 0.0061 ## collgrad 0.68770 0.0077 ## hours 0.13666 0.0059 ## ##$percentual contribution per dimension
##          dim. % contribution     SE
## wage                 0.32908 0.0083
## $percentual contribution per dimension ## dim. % contribution SE ## wage 0.33333 0 ## collgrad 0.33333 0 ## hours 0.33333 0 # PDF page 22 for (ks in seq(1 / 3, 1, 1 / 3)) { print( svyafcdec( ~ wage + collgrad + hours, subgroup = ~ factor(married), design = des_nlsw88, cutoffs = list(4, 'college grad' , 26), k = ks , g = 0, dimw = NULL, na.rm = T )[6] ) } ##$percentual contribution per subgroup
## $percentual contribution per subgroup ## grp. % contribution SE ## 0 0.32032 0.0197 ## 1 0.67968 0.0197 ## ##$percentual contribution per subgroup
## 1              0.6701 0.0477