Module C2 Track B lab (Part 2 of 3) Distribution: extend the open age group and smooth
Reference solution (full code) for expert review
13 August 2026
The lab census is single years of age topped at 75+. A base population has to reach 100+, and its single-year counts have to be believable — free of the spikes people’s rounded ages put there. This part does both: it extends the open age group to 105+, and it decides, by the protocol’s published rule, how strong a moving average the heaping deserves. Those are steps D and E of the Population Division’s method protocol.
Because the series is already single years, nothing is graduated here: the protocol’s single-year route measures the heaping directly and smooths the counts as they stand.
You do not need to have done Part 1. Sheet 01 of this workbook is the table Part 1 produces — the census counts already corrected for coverage — so this part starts from a known place. Fill in rows D and E of the decision log (sheet 04) as you go.
Running it. Put UGA_C2_lab2_data.xlsx
in the same folder as this file. Then either run the chunks one at a
time in RStudio, or press Knit to produce this report —
from the console that is
rmarkdown::render("C2_lab2_distribution.Rmd"). Nothing else
is needed: the protocol functions travel with the file, in an appendix
at the end.
library(readxl) # reads the workbook
library(dplyr) # data handling the workflow functions use
library(tibble) # tibble(), used by one of the protocol functions
library(DemoTools) # the demographic methods underneath everything
library(knitr) # kable(), for readable tables
DATA <- "UGA_C2_lab2_data.xlsx"
## Stop early, with directions, if the workbook cannot be opened: a missing
## file, the wrong working directory, or a OneDrive file that has not been
## downloaded to this computer all land here.
if (!isTRUE(tryCatch({ excel_sheets(DATA); TRUE }, error = function(e) FALSE)))
stop("Cannot open ", DATA, " from ", getwd(),
". Put the workbook in the same folder as this file, or give DATA ",
"its full path. If it sits in OneDrive, right-click the file, choose ",
"'Always keep on this device', close it in Excel, and try again.")1 Read the workbook
start <- read_excel(DATA, sheet = "01_start", skip = 6, n_max = 76)
lt <- read_excel(DATA, sheet = "02_lifetable", skip = 3, n_max = 19)
cov <- read_excel(DATA, sheet = "03_covariates", skip = 2, n_max = 1)
details <- read_excel(DATA, sheet = "01_start", range = "A4:B5",
col_names = c("field", "value"))
REFDATE <- as.numeric(details$value[details$field == "census_reference_date"])
LOCID <- as.numeric(details$value[details$field == "LocID"])One setting that is a choice rather than data. It belongs in the decision log, not in the workbook, which is why it sits here in plain sight.
OANEW <- 105 # extend to this open age group2 Step D — extend the open age group
75+ is a box with no shape inside it — three-quarters of a million people with no ages. The extension replaces it with single years running to 105+, using the life table as the pattern: it fits a stable population to the ages below the box, then asks how the people inside it would be spread if they died at the life table’s rates. The function also looks for the age above which the reported counts stop behaving — old-age exaggeration usually reaches below the open group — and redistributes from there up, so it repairs more than the box alone.
ext <- census_workflow_extend_OAG(popM = start$Male, popF = start$Female, Age = start$AgeStart,
LocID = LOCID, Year = floor(REFDATE),
lxM = lt$lx_M, lxF = lt$lx_F, Age_lx = lt$Age,
OAnew = OANEW)
extended <- data.frame(Age = ext$Age_ext, Male = ext$popM_ext, Female = ext$popF_ext)
cat("redistribution starts at age", ext$age_redist_start,
"| open age group is now", max(extended$Age), "\n")## redistribution starts at age 65 | open age group is now 105pstart()
old <- start[start$AgeStart >= 40 & start$AgeStart < 75, ]
new <- extended[extended$Age >= 40, ]
plot(new$Age, (new$Male + new$Female) / 1000, type = "l", lwd = 2, col = BLUE, bty = "l",
xlab = "Age", ylab = "Thousands per single year of age", xlim = c(40, 105))
lines(old$AgeStart, (old$Male + old$Female) / 1000, lwd = 2, col = GREY)
abline(v = ext$age_redist_start, lty = 3, col = GREY)
text(ext$age_redist_start, par("usr")[4], paste(" redistribution from", ext$age_redist_start),
adj = c(0, 1.4), col = GREY, cex = 0.9)
legend("topright", c("Reported single years (to 74)", "Extended to 105+"),
lwd = 2, col = c(GREY, BLUE), bty = "n")
Note where the redistribution starts: below the open group, not at it. The counts between that age and 75 were reported, but the extension replaces them too, because age exaggeration at the oldest ages reaches down past the box. Everything below the start age is untouched.
Your turn. Check the arithmetic: sum the extended series from the redistribution start age upward and compare it with the same range in the series you started from. They should match — the extension redistributes, it does not add. Record the start age in the decision log.
3 Step E — smooth the single years
No graduation happens here: the series is already single years. What remains is the heaping Part 1 measured. The protocol does not smooth by taste: it reads the Bachi index, looks at how concentrated the heaping is on the favourite digits, reads mean years of schooling as a proxy for how well ages are reported, and picks a moving-average width from a published rule — separately for children and for adults, because the two are disturbed by different things.
sm <- census_workflow_adjust_smooth(popM = extended$Male, popF = extended$Female,
Age = extended$Age,
EduYrs = cov$Value[1])
kable(data.frame(
Series = c("Adults", "Children"),
`Bachi index` = round(c(sm$bachi_adult, sm$bachi_child), 2),
`Mean years of schooling` = cov$Value[1],
`Smoothing the rule chose` = c(sm$best_smooth_adult, sm$best_smooth_child),
check.names = FALSE))| Series | Bachi index | Mean years of schooling | Smoothing the rule chose |
|---|---|---|---|
| Adults | 3.85 | 6.3 | bestMavN = 4 |
| Children | 1.81 | 6.3 | bestMavN = 2 |
Read the table across. bestMavN = 4 means a four-year
moving average — wide enough to absorb spikes two years either side of a
round age. The children’s window is narrower because their ages,
reported by parents of recently born children, carry less rounding.
pstart()
plot(0:105, (sm$popM_smoothed + sm$popF_smoothed) / 1000, type = "l", lwd = 2, col = BLUE,
bty = "l", xlab = "Age", ylab = "Thousands per single year of age")
lines(extended$Age, (extended$Male + extended$Female) / 1000, lwd = 1.2, col = GREY)
legend("topright", c("Smoothed", "Before smoothing"), lwd = c(2, 1.2),
col = c(BLUE, GREY), bty = "n")
cat("total before", format(round(sum(extended$Male + extended$Female)), big.mark = ","),
"| after", format(round(sum(sm$popM_smoothed + sm$popF_smoothed)), big.mark = ","), "\n")## total before 47,968,040 | after 47,968,040The two totals match because the function rescales the smoothed series back to the totals it was given — the protocol’s proration step. Smoothing moves people between ages; it must never create or lose any.
Your turn. Run this part as it stands, then change EduYrs on sheet 03 from 6.3 to 9 and run it again. The rule narrows both windows — better age reporting deserves lighter smoothing — and the printed table shows it. Record both outcomes in the decision log.
dir.create("out2", showWarnings = FALSE)
write.csv(data.frame(Age = 0:105, Male = sm$popM_smoothed, Female = sm$popF_smoothed),
"out2/part2_smoothed_single_year.csv", row.names = FALSE)
write.csv(extended, "out2/part2_unsmoothed_single_year.csv", row.names = FALSE)
writeLines(sm$best_smooth_child, "out2/part2_smooth_method.txt")4 Where this part ends
The population is now single years of age running to 105+, smoothed, at the level the coverage adjustment set. What it still carries is the census’s shortage of young children, which no amount of smoothing can repair, because the children were never counted. Part 3 rebuilds them from births and survivorship, and dates the result.
5 Bring your own census
Copy the provided input file, clear the green cells, and enter your own values with their sources and reference dates.
Expect that to take much longer than the lab does, and to happen outside the session. Many offices may find that some of the inputs are not published for their country, and finding out which ones is worth as much as the run itself. Record each one in the decision log.
6 Appendix — the protocol functions this part uses
These are the Population Division’s own censusAdjust
files.Underneath they call DemoTools; what they add is the protocol’s
decision rules.
## ==== census_workflow_extend_OAG.R ==========================
# Extend and redistribute open age group. Default used for R22 is 105+
census_workflow_extend_OAG <- function(popM,
popF,
Age,
LocID,
Year,
cv_tolerance = 0.75,
min_age_redist = min(max(Age), 65),
OAnew = 105,
lxM = NULL,
lxF = NULL,
Age_lx = NULL,
AgeInt_lx = NULL) {
maxage <- max(Age)
stopifnot(maxage >= 50 & OAnew > maxage)
single <- DemoTools::is_single(Age)
# group to five-year age groups
popM5 <- DemoTools::groupAges(popM, Age=Age, N=5)
popF5 <- DemoTools::groupAges(popF, Age=Age, N=5)
Age5 <- seq(0, maxage, 5)
# If OAnew is not a multiple of 5, round up
OpenAge5 <- as.integer(OAnew / 5) * 5
OAnew5 <- ifelse(OAnew == OpenAge5, OAnew, OpenAge5+5)
# extend pop to OAnew using OPAG
if(popM5[Age5 == maxage] > 5) { # male pop at older ages tends to be smaller than female, so use this as criteria for extension
if (is.null(lxM)) {
# download abridged life table for males
ltM_abr <- fertestr::FetchLifeTableWpp2019(locations = LocID,
year = Year,
sex = "male")
# download abridged life table for males
ltF_abr <- fertestr::FetchLifeTableWpp2019(locations = LocID,
year = Year,
sex = "female")
lxM <- ltM_abr$lx
lxF <- ltF_abr$lx
Age_lx <- ltM_abr$x
}
# keep only lx for abridged age groups
ages_abr <- c(0,1,seq(5, max(Age_lx), 5))
lxM <- lxM[Age_lx %in% ages_abr]
lxF <- lxF[Age_lx %in% ages_abr]
# recompute life table, extending as necessary
nLxM <- lt_abridged(lx = lxM, Age = ages_abr, Sex = "m", OAnew = OAnew5)$nLx
nLxF <- lt_abridged(lx = lxF, Age = ages_abr, Sex = "f", OAnew = OAnew5)$nLx
# collapse first two age groups
nLxM <- c(nLxM[1]+nLxM[2],nLxM[3:length(nLxM)])
nLxF <- c(nLxF[1]+nLxF[2],nLxF[3:length(nLxF)])
Age_nLx <- seq(0,OAnew5,5)
redist_ages <- seq(min_age_redist, maxage, 5)
nra <- length(redist_ages)
pop_extM <- list()
pop_extF <- list()
cvM <- list()
cvF <- list()
# extend the old-age populations using different starting ages
for (k in 1:nra) {
pop_extM[[k]] <- OPAG(Pop = popM5,
Age_Pop = Age5,
nLx = nLxM,
Age_nLx = Age_nLx,
method = "mono",
Redistribute_from = redist_ages[k],
OAnew = OAnew5)
cvM[[k]] <- sd(diff(pop_extM[[k]]$Pop_out[seq(0,OAnew5,5) %in% c((redist_ages[k]-10):(redist_ages[k]+5))]))/abs(mean(diff(pop_extM[[k]]$Pop_out[seq(0,OAnew5,5) %in% c((redist_ages[k]-10):(redist_ages[k]+5))])))
pop_extF[[k]] <- OPAG(Pop = popF5,
Age_Pop = Age5,
nLx = nLxF,
Age_nLx = Age_nLx,
method = "mono",
Redistribute_from = redist_ages[k],
OAnew = OAnew5)
cvF[[k]] <- sd(diff(pop_extF[[k]]$Pop_out[seq(0,OAnew5,5) %in% c((redist_ages[k]-10):(redist_ages[k]+5))]))/abs(mean(diff(pop_extF[[k]]$Pop_out[seq(0,OAnew5,5) %in% c((redist_ages[k]-10):(redist_ages[k]+5))])))
}
cvM <- do.call(rbind, cvM)
cvF <- do.call(rbind, cvF)
cvM_below_tolerance <- min(cvM) <= cv_tolerance
cvF_below_tolerance <- min(cvF) <= cv_tolerance
if (cvM_below_tolerance & cvF_below_tolerance) {
# identify the max age at which the cv tolerance is below the threshold for both males and females
age_redist <- min(max(redist_ages[cvM <= cv_tolerance]), max(redist_ages[cvF <= cv_tolerance]))
} else {
# use the minimum of the age with minimum cv for males and age with minimum cv for females
age_redist <- min(redist_ages[cvM == min(cvM)], redist_ages[cvF == min(cvF)])
}
Agein <- Age
Age <- pop_extM[[1]]$Age_out
popM_ext <- pop_extM[[c(1:nra)[redist_ages==age_redist]]]$Pop_out
popF_ext <- pop_extF[[c(1:nra)[redist_ages==age_redist]]]$Pop_out
############
############
############
} else { # if male pop in maxage group is less than 5, then extend both males and females to to OAnew5 with zeros
age_redist <- maxage
cv_out <- NA
popM_ext <- c(popM5, rep(0, (OAnew5 - maxage) / 5))
popF_ext <- c(popM5, rep(0, (OAnew5 - maxage) / 5))
Agein <- Age
Age <- seq(0,OAnew5,5)
}
if (single) { # if input is single, then graduate output
# Now graduate to single year of age
popM_ext <- DemoTools::graduate_mono(popM_ext, AgeInt = rep(5,(OAnew5 + 5)/5), Age= seq(0,OAnew5,5), OAG = TRUE)
popF_ext <- DemoTools::graduate_mono(popF_ext, AgeInt = rep(5,(OAnew5 + 5)/5), Age= seq(0,OAnew5,5), OAG = TRUE)
Age <- 0:OAnew5
# truncate back to the originally requested OAnew
popM_ext <- DemoTools::groupOAG(Value = popM_ext, Age = Age, OAnew = OAnew)
popF_ext <- DemoTools::groupOAG(Value = popF_ext, Age = Age, OAnew = OAnew)
Age <- 0:OAnew
}
# splice together, using opag extension for only ages above the start of the opag redistribution
popM_ext_final <- c(popM[Agein < age_redist],
popM_ext[Age >= age_redist])
popF_ext_final <- c(popF[Agein < age_redist],
popF_ext[Age >= age_redist])
Age_ext_final <- c(Agein[Agein < age_redist],
Age[Age >= age_redist])
opag_out <- list(popM_ext = popM_ext_final,
popF_ext = popF_ext_final,
Age_ext = Age_ext_final,
age_redist_start = age_redist)
return(opag_out)
}
## ==== census_workflow_adjust_smooth.R =======================
## census_workflow_adjust_smooth
##
## this function assesses age heaping and accordingly applies smoothing
## as described in protocol flowchart xx
##
## @param popM numeric. vector of male population counts
## @param popF numeric. vector of female population counts
## @param Age numeric. vector of age (years) associated with population counts vectors
## @return A list of objects describing the steps implemented in smoothing population counts over age
## @export
census_workflow_adjust_smooth <- function(popM,
popF,
Age,
bachi_age_child = 3:17, # age range for bachi index for children
bachi_age_adult = 23:77, # age range for bachi index for adults
age_ratio_age_child = c(0,10), # age range for age ratio score for children
age_ratio_age_adult = c(15,70), # age range for age ratio score for adults
EduYrs, # average years of education (used as a criterion for level of smoothing)
graduation_method = "graduate_mono") {
# intialize bachi
bachi_child <- NA
bachi_adult <- NA
# if inputs are by single year of age
if (is_single(Age)) {
# assess single year age heaping for children and smooth accordingly
pop_smooth_child <- getSmoothedPop1 (Age = Age,
popF = popF,
popM = popM,
bachi_age = bachi_age_child,
age_ratio_age = age_ratio_age_child,
EduYrs = EduYrs,
subgroup = "child")
bachi_child <- pop_smooth_child$bachi
# assess single year age heaping for adults and smooth accordingly
pop_smooth_adult <- getSmoothedPop1 (Age = Age,
popF = popF,
popM = popM,
bachi_age = bachi_age_adult,
age_ratio_age = age_ratio_age_adult,
EduYrs = EduYrs,
subgroup = "adult")
bachi_adult <- pop_smooth_adult$bachi
} else {
# assess grouped age heaping for children and smooth accordingly
pop_smooth_child <- getSmoothedPop5(Age = Age,
popF = popF,
popM = popM,
age_ratio_age = age_ratio_age_child,
EduYrs = EduYrs,
subgroup = "child",
graduation_method = graduation_method)
# assess grouped age heaping for adults and smooth accordingly
pop_smooth_adult <- getSmoothedPop5(Age = Age,
popF = popF,
popM = popM,
age_ratio_age = age_ratio_age_adult,
EduYrs = EduYrs,
subgroup = "adult",
graduation_method = graduation_method)
}
# blend the smoothed child and adult series, with transition at ages 15-19
wts <- c(rep(1,16),0.8, 0.6, 0.4, 0.2, rep(0, max(Age)-19))
popM_smoothed <- (pop_smooth_child$popM_smooth * wts) + (pop_smooth_adult$popM_smooth * (1-wts))
popF_smoothed <- (pop_smooth_child$popF_smooth * wts) + (pop_smooth_adult$popF_smooth * (1-wts))
# re-adjust to ensure that after smoothing we are still matching the total
popM_smoothed <- popM_smoothed * sum(popM)/sum(popM_smoothed)
popF_smoothed <- popF_smoothed * sum(popF)/sum(popF_smoothed)
pop_smoothed <- list(Age = 0:105,
popF_smoothed = popF_smoothed,
popM_smoothed = popM_smoothed,
bachi_child = bachi_child,
bachi_adult = bachi_adult,
ageRatio_adult_orig = pop_smooth_adult$AgeRatioScore_orig,
ageRatio_child_orig = pop_smooth_child$AgeRatioScore_orig,
ageRatio_adult_mav2 = pop_smooth_adult$AgeRatioScore_mav2,
ageRatio_child_mav2 = pop_smooth_child$AgeRatioScore_mav2,
best_smooth_adult = pop_smooth_adult$best_smooth_method,
best_smooth_child = pop_smooth_child$best_smooth_method)
return(pop_smoothed)
}
## ==== census_workflow_getSmoothedPop1.R =====================
# for single year data, assesses heaping with bachi index and smoothes with moving avearage at different levels
# for abridged or five-year data, smoothes with moving average at different levels
# inputs are population by single year of age and sex
getSmoothedPop1 <- function(popM,
popF,
Age,
bachi_age = NULL,
age_ratio_age = NULL, # must be multiples of 5
EduYrs,
subgroup = c("adult", "child")) {
maxage = max(Age)
if (is.null(bachi_age)) {
if (subgroup == "adult") {
ageMin <- 23
ageMax <- min(77, maxage)
} else {
ageMin <- 3
ageMax <- 17
}
} else {
ageMin <- min(bachi_age)
ageMax <- max(bachi_age)
}
# compute bachi
bachi_m <- check_heaping_bachi(Value = popM,
Age = Age,
ageMin = ageMin,
ageMax = ageMax,
method = "pasex",
details = TRUE)
bachi_f <- check_heaping_bachi(Value = popF,
Age = Age,
ageMin = ageMin,
ageMax = ageMax,
method = "pasex",
details = TRUE)
# compute proportion of heaping concentrated in digits 0 and 5
BachiProp0and5_m <- (bachi_m$pct[1] + bachi_m$pct[6] - 20) / bachi_m$index
BachiProp0and5_f <- (bachi_f$pct[1] + bachi_f$pct[6] - 20) / bachi_f$index
# compute proportion of heaping concentrated in the favorite two digits
BachiPropMax2_m <- (sum(bachi_m$pct[order(-bachi_m$pct)][c(1,2)]) - 20) / bachi_m$index
BachiPropMax2_f <- (sum(bachi_f$pct[order(-bachi_f$pct)][c(1,2)]) - 20) / bachi_f$index
# identify the preferred smoothing based on bachi, digit preference and level of education
bestMavN <- getBestMavN(Bachi = max(bachi_m$index, bachi_f$index),
BachiProp0and5 = min(BachiProp0and5_m, BachiProp0and5_f),
BachiPropMax2 = min(BachiPropMax2_m, BachiPropMax2_f),
EduYrs = EduYrs,
subgroup = subgroup)
# smooth based on bestMavN
if (!is.na(bestMavN)) {
popM_smooth_mav <- mavPop1(popM, Age = Age)
popF_smooth_mav <- mavPop1(popF, Age = Age)
# Parse the best single age data
popM_smooth <- unlist(select(popM_smooth_mav$MavPopDF, !!paste0("Pop", bestMavN)))
popF_smooth <- unlist(select(popF_smooth_mav$MavPopDF, !!paste0("Pop", bestMavN)))
best_smooth_method <- paste("bestMavN = ", bestMavN)
AgeRatioScore_orig <- NA # We dont need age ratio score if we are using single age data so set this to NULL for output
AgeRatioScore_mav2 <- NA
} else { # if bestMavN is NA then group to five year data
# group single year series to five year age groups
popM5 <- DemoTools::groupAges(popM, Age = Age, N = 5, OAnew = maxage)
popF5 <- DemoTools::groupAges(popF, Age = Age, N = 5, OAnew = maxage)
Age5 <- seq(0, maxage, 5)
nAge5 <- length(Age5)
# smooth the 5-year series using mav2
popM5_mav2 <- DemoTools::smooth_age_5(popM5, seq(0, maxage, 5), method = "MAV", n = 2)
popF5_mav2 <- DemoTools::smooth_age_5(popF5, seq(0, maxage, 5), method = "MAV", n = 2)
# smooth the 5-year series using mav4
popM5_mav4 <- DemoTools::smooth_age_5(popM5, seq(0, maxage, 5), method = "MAV", n = 4)
popM5_mav4[2] <- popM5_mav2[2]
popM5_mav4[nAge5 - 2] <- popM5_mav2[nAge5 - 2]
popF5_mav4 <- DemoTools::smooth_age_5(popF5, seq(0, maxage, 5), method = "MAV", n = 4)
popF5_mav4[2] <- popF5_mav2[2]
popF5_mav4[nAge5 - 2] <- popF5_mav2[nAge5 - 2]
# compute age ratio scores for 5-yr age groups from 15-19 to 70-74 for adults or from 0 to 10-14 for children
if (is.null(age_ratio_age)) {
if (subgroup == "adult") {
ageMin <- 15
ageMax <- min(70, maxage)
} else {
ageMin <- 0
ageMax <- 10
}
} else {
ageMin = min(age_ratio_age)
ageMax = max(age_ratio_age)
}
# first on the unsmoothed 5-year data
ageRatioScoreM_orig <- DemoTools::ageRatioScore(Value = popM5, Age = Age5,
ageMin = ageMin, ageMax = ageMax, OAG = FALSE)
ageRatioScoreF_orig <- DemoTools::ageRatioScore(Value = popF5, Age = Age5,
ageMin = ageMin, ageMax = ageMax, OAG = FALSE)
AgeRatioScore_orig = max(ageRatioScoreM_orig, ageRatioScoreF_orig)
# then on the mav2 smoothed 5-year data
ageRatioScoreM_mav2 <- DemoTools::ageRatioScore(Value = popM5_mav2, Age = Age5,
ageMin = ageMin, ageMax = ageMax, OAG = FALSE)
ageRatioScoreF_mav2 <- DemoTools::ageRatioScore(Value = popF5_mav2, Age = Age5,
ageMin = ageMin, ageMax = ageMax, OAG = FALSE)
AgeRatioScore_mav2 = max(ageRatioScoreM_mav2, ageRatioScoreF_mav2)
# now identify the best smoothing approach based on age ratio scores and education
bestGrad5 <- getBestGrad5(AgeRatioScore_orig = AgeRatioScore_orig,
AgeRatioScore_mav2 = AgeRatioScore_mav2,
EduYrs = EduYrs,
subgroup = subgroup)
if (bestGrad5 == 1) {
popM_smooth <- DemoTools::graduate_mono(popM5, AgeInt = DemoTools::age2int(Age5), Age = Age5, OAG = TRUE)
popF_smooth <- DemoTools::graduate_mono(popF5, AgeInt = DemoTools::age2int(Age5), Age = Age5, OAG = TRUE)
}
if (bestGrad5 == 2) {
popM_smooth <- DemoTools::graduate_mono(popM5_mav2, AgeInt = DemoTools::age2int(Age5), Age = Age5, OAG = TRUE)
popF_smooth <- DemoTools::graduate_mono(popF5_mav2, AgeInt = DemoTools::age2int(Age5), Age = Age5, OAG = TRUE)
}
if (bestGrad5 == 4) {
popM_smooth <- DemoTools::graduate_mono(popM5_mav4, AgeInt = DemoTools::age2int(Age5), Age = Age5, OAG = TRUE)
popF_smooth <- DemoTools::graduate_mono(popF5_mav4, AgeInt = DemoTools::age2int(Age5), Age = Age5, OAG = TRUE)
}
best_smooth_method <- paste("bestGrad5 = ", bestGrad5)
} # done with smoothing
out.data <- list(popM_smooth = popM_smooth,
popF_smooth = popF_smooth,
best_smooth_method = best_smooth_method,
bachi = max(bachi_m$index, bachi_f$index),
AgeRatioScore_orig = AgeRatioScore_orig,
AgeRatioScore_mav2 = AgeRatioScore_mav2)
return(out.data)
}
## ==== census_workflow_getBestMavN.R =========================
#Get BestMavN from Bachi Indices ---------------------------------------------
getBestMavN <- function(Bachi, BachiProp0and5, BachiPropMax2, EduYrs, subgroup = c("adult", "child")) {
## NEED TO IMPLEMENT DIFFERENT FILTERING BASED ON SUBGROUP????
# select best MavN
BachiLevels <- c(0, 0.75, 2, 4, 8, 30, Inf)
MavNs <- c( 1, 2, 4, 6, 10, NA)
InitialBest <- MavNs[as.double(cut(Bachi, breaks = BachiLevels, labels = FALSE))]
FinalBest <- InitialBest
if (!is.na(FinalBest)) {
# these could be nested if else.. but this is easier to read
if (InitialBest == 6) {
if(BachiProp0and5 > 0.65) {
FinalBest <- 10} else {
if (EduYrs >= 8) FinalBest <- 4
if (EduYrs < 8) FinalBest <- 6
}
}
if (InitialBest == 4) {
if(BachiProp0and5 > 0.60) {
FinalBest <- 6 } else {
if (EduYrs >= 8) FinalBest <- 2
if (EduYrs < 8) FinalBest <- 4
}
}
if (InitialBest == 2) {
if(BachiPropMax2 > 0.70) {
FinalBest <- 4} else {
if (EduYrs >= 8) FinalBest <- 1
}
}
if (InitialBest == 1) {
if(BachiPropMax2 > 0.55) {
FinalBest <- 2}
}
} # close for if !is.na(FinalBest)
return(FinalBest)
}
## ==== census_workflow_getBestGrad5.R ========================
#Get BestGrad5 from smoothing results
getBestGrad5 <- function(AgeRatioScore_orig, AgeRatioScore_mav2, EduYrs, subgroup = c("adult","child")) {
if (subgroup == "adult") {
# select whether to use the straight 5-year data
# or use Mav2 or mav4 of the 5-year data
BestGrad5 <- NA
if (AgeRatioScore_orig < 4 ) {
BestGrad5 <- 1
} else {
if (EduYrs >= 4) {
BestGrad5 <-1
} else {
if (AgeRatioScore_mav2 < 4) {
BestGrad5 <- 2
} else {
BestGrad5 <- 4
}
}
}
}
if (subgroup == "child") {
BestGrad5 <- NA
if (AgeRatioScore_orig < 4 ) {
BestGrad5 <- 1
} else {
if (EduYrs >= 4) {
BestGrad5 <-1
} else {
BestGrad5 <- 2
}
}
}
return(BestGrad5)
}
## ==== getmavPop1.R ==========================================
# Test of functions to create moving average estimates of population by single years of age
# and birth cohort estimates based on single year data and reference date
# note Bachi defaults ageMin=23 and ageMax=77
# 11/19/20 rc code to fill in youngest oldest ages using seq of mavs
# 11/25/20 fixed mutate
# 12/17/20 more fixes
# 2021-01-20 remove BirthCohort1 (switch to DT function )
# No longer do birth cohorts here
# 2021-01-29 remove sex, CenDate, and CenYr entries from mavPop1
#
SCcount <- function(Pop) {
# count number of sign changes age to age
# get Pop(a) - Pop(a-1) only look at ages 23-67
PopCh <- Pop[25:68] - Pop[24:67]
PopCHsign <- sign(PopCh)
PopCHsignChange <- PopCHsign[1:43] != PopCHsign[2:44]
sum(PopCHsignChange)
}
#
# new version of function with help from rc 11/19/20
mavPop1 <- function(Pop, Age) {
# store pop data for sending back
MavPopDF <- tibble(Age = Age, Pop1 = Pop)
Nage <- length(Age) - 2
# get results for original data
BachiN <- check_heaping_bachi(Pop, Age, method = "pasex")
MySCcount <- SCcount(Pop)
MavSummaryDF <- tibble(MavN = 1,BachiN, SCcountN = MySCcount)
for (MavN in 2:10) {
MyMavNew <- mav(Pop, n = MavN, Age = Age)
MavPopDF <- cbind(MavPopDF, MyMavNew)
colnames(MavPopDF)[ncol(MavPopDF)] <- paste("Pop", MavN, sep = "")
}
OpenAge <- Nage + 1
MavPopDF <- mutate( MavPopDF, across(3:11, ~(
function(x) {
for(i in 1:length(x))
{
if(is.na(x[i]) & (Age[i] == 0 | Age[i] >= (OpenAge - 1)))
x[i] <- Pop1[i]
else if(is.na(x[i]) & (Age[i] == 1 | Age[i] == (OpenAge - 2)))
x[i] <- Pop2[i]
else if(is.na(x[i]) & (Age[i] == 2 | Age[i] == (OpenAge - 3)))
x[i] <- Pop4[i]
else if(is.na(x[i]) & (Age[i] == 3 | Age[i] == (OpenAge - 4)))
x[i] <- Pop6[i]
else if(is.na(x[i]) & (Age[i] == 4 | Age[i] == (OpenAge - 5)))
x[i] <- Pop8[i]
}
return(x)
} ) # end of function spec
(.) ) # end of across spec I still don't know the function of (.)
) # end of mutate
for (i in 3:ncol(MavPopDF)){
BachiN <- check_heaping_bachi(MavPopDF[[i]], Age, method = "pasex")
SCcountN <- SCcount(MavPopDF[[i]])
MavSummaryDF <- rbind(MavSummaryDF,c(MavN=i, BachiN, SCcountN))
}
rownames(MavPopDF) <- MavPopDF$Age
list(MavPopDF = MavPopDF,
MavSummaryDF = MavSummaryDF)
} # end of mavPop1 function