Análisis demográfico en R: técnicas básicas y DemoTools

Introducción

Objetivo y setting inicial

El objetivo de este espacio es recorrer algunas técnicas básicas para el análisis demográfico, específicamente en lo que nos toca como profesionales al levantar cada censo. Utlizaremos funciones ya implementadas en paquetes de R y otras que no, pudiendo los asistentes adaptarlas a sus propios fines con posterioridad al workshop. Trataremos ejemplos prácticos con información de países de América Latina.

Nos valdremos principalmente del paquete DemoTools, elaborado por Tim Riffe y colaboradores. ¡Gracias!

¿Qué formato tendremos? Lo pensamos como taller, por lo que iremos corriendo el código a la par, tratando de que nos sigamos mutuamente. Si queda algo picando, podemos seguir conectándonos a través de correo electrónico.

Daremos por sentado cierto conocimiento de R y específicamente de dplyr and ggplot. Si no sos del caso te recomendamos ver algunos de los siguientes materiales: web del paquete, o una viñeta introductoria.

Se requiere la instalación de los siguientes paquetes:

install.packages("tidyverse")
install.packages("readxl")
install.packages("mipfp")
install.packages("LexisPlotR")
install.packages("MortalityLaws")
install.packages("remotes")
install.packages("quadprog")
install.packages("kableExtra")
install.packages("MortCast")
install.packages("zoo")

# y la instalación de DemoTools requiere dos pasos
install.packages("rstan", repos = c("https://mc-stan.org/r-packages/", getOption("repos")))
remotes::install_github("timriffe/DemoTools")

El repositorio donde esta alojado esta web lo podes encontrar aquí. Esperamos que este contenido les sirva para complementar su trabajo en el contexto de esta ronda censal.

Data

Como primer paso, creemos un proyecto R con todo el contenido descargado en el mismo directorio. Activemos las librerías:

library(tidyverse)
library(DemoTools)
library(kableExtra)

En este enlace pueden encontrar un archivo data_LA.rda con información sobre la población por edad y sexo relevada en censos de países de Latinoaméríca y el Caribe, defunciones y nacimientos por edad y sexo, recolectadas las tres series a partir de las publicaciones Demographic Yearbook de Naciones Unidas, y disponibles en el portal de datos DemoData. Demos un primer vistazo a la información:

load("data_LA.rda")
glimpse(census_LA)

## Rows: 25,827
## Columns: 11
## $ LocID       <int> 32, 32, 32, 32, 32, 32, 32, 32, 32, 32, 32, 32, 32, 32, 32…
## $ País        <chr> "Argentina", "Argentina", "Argentina", "Argentina", "Argen…
## $ Año         <int> 1960, 1960, 1960, 1960, 1960, 1960, 1960, 1960, 1960, 1960…
## $ Año_medio   <dbl> 1960.75, 1960.75, 1960.75, 1960.75, 1960.75, 1960.75, 1960…
## $ Tipo        <chr> "De-facto", "De-facto", "De-facto", "De-facto", "De-facto"…
## $ Fecha       <chr> "30/09/1960", "30/09/1960", "30/09/1960", "30/09/1960", "3…
## $ Sexo        <chr> "Hombre", "Hombre", "Hombre", "Hombre", "Hombre", "Hombre"…
## $ Edad        <int> 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, …
## $ Edad_interv <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
## $ Edad_grupo  <chr> "0", "1", "2", "3", "4", "5", "6", "7", "8", "9", "10", "1…
## $ Pop         <dbl> 217996, 213352, 214570, 217661, 215289, 212812, 209918, 21…

La cobertura geogrpafica y temporal de la inforamción censal inlcuye inforamción sobre censos recientes (Panamá, Brasil y Ecuador):

census_LA %>% 
     distinct(Año, País) %>% 
     ggplot(aes(Año, País)) + 
     geom_point() +
     ggtitle("Población por edad y sexo de censos \nincluidos para el workshop") +
     theme_minimal()

Respecto a las defunciones y nacimientos, las tablas deaths_LA y births_LA contienen información de eventos por sexo (del recién nacido en el primer caso) y edad (de madres en el primer caso). Por ejemplo Panamá en 2022:

deaths_LA %>% 
     filter(País == "Panama", Año == 2022) %>% 
     head()

##   LocID   País  Año Año_medio               Tipo Edad   Sexo Edad_interv
## 1   591 Panama 2022    2022.5 Year of occurrence    0 Hombre           1
## 2   591 Panama 2022    2022.5 Year of occurrence    1 Hombre           4
## 3   591 Panama 2022    2022.5 Year of occurrence    5 Hombre           5
## 4   591 Panama 2022    2022.5 Year of occurrence   10 Hombre           5
## 5   591 Panama 2022    2022.5 Year of occurrence   15 Hombre           5
## 6   591 Panama 2022    2022.5 Year of occurrence   20 Hombre           5
##   Edad_grupo Deaths
## 1       <NA>    450
## 2       <NA>    134
## 3       <NA>     40
## 4       <NA>     58
## 5       <NA>    161
## 6       <NA>    247

births_LA %>% 
     filter(País == "Panama", Año == 2022) %>% 
     head()

##   LocID   País  Año Año_medio               Tipo   Sexo Edad Edad_interv
## 1   591 Panama 2022    2022.5 Year of occurrence Hombre   10           5
## 2   591 Panama 2022    2022.5 Year of occurrence Hombre   15           5
## 3   591 Panama 2022    2022.5 Year of occurrence Hombre   20           5
## 4   591 Panama 2022    2022.5 Year of occurrence Hombre   25           5
## 5   591 Panama 2022    2022.5 Year of occurrence Hombre   30           5
## 6   591 Panama 2022    2022.5 Year of occurrence Hombre   35           5
##   Edad_grupo Births
## 1       <NA>    206
## 2       <NA>   4653
## 3       <NA>   9030
## 4       <NA>   8154
## 5       <NA>   6130
## 6       <NA>   3489

bind_rows(deaths_LA %>% mutate(Evento = "Defunciones"),
          births_LA %>% mutate(Evento = "Nacimientos")) %>% 
     distinct(Año, País, Evento) %>% 
     ggplot(aes(Año, País)) + 
     geom_line(data = . %>% filter(Evento == "Nacimientos")) +
     geom_point(data = . %>% filter(Evento == "Defunciones"), col = 2) +
     ggtitle("Nacimientos (línea) y defunciones (punto rojo) \nincluidos para el workshop") +
     theme_minimal()

Si quieres ver en profundidad la obtención de la información puedes ver el script get_data.R en el repositorio.

Utilizaremos funciones que no estan implementadas aún en paquetes, incluidas en el script:

source("funs.R")

Población

Evolución censal

Sigamos con el caso de Panamá, aprovechando la publicación reciente sobre el censo 2023. Fue un censo de-jure, con campo del 8 de enero al 6 de marzo del 2023, resultando en más de 4 millones de personas residentes. La dinámica demográfica de este país estuvo muy relacionada con el desarrollo del canal, influenciando el nivel de desarrollo económico y la condición de país atractor para la migración (Estadística y Censos 2016).

Resumamos en la siguiente tabla la evolución de la población censada, la fecha y tipo de censo, el ritmo de crecimiento (tasa de crecimiento medio anual), el porcentaje de mujeres y la edad media:

panama <- census_LA %>% filter(País == "Panama")
panama %>% 
     summarise(Poblacion = sum(Pop),
               Tipo = unique(Tipo),
               Fecha = dmy(unique(Fecha)),
               Porc_Muj = round(sum(Pop[Sexo == "Mujer"])/Poblacion, 3),
               Edad_media = round(sum(Edad * Pop, na.rm = T) / sum(Pop), 3),
               .by = Año) %>% 
     mutate(tca = round(log(lead(Poblacion)/Poblacion)/(lead(dec.date(Fecha))-dec.date(Fecha)), 3)) %>% 
     kbl(caption = "Evolución de indicadores censales. Panamá") %>% 
     kable_classic(full_width = F)

Evolución de indicadores censales. Panamá
Año	Poblacion	Tipo	Fecha	Porc_Muj	Edad_media	tca
1950	758562	De-facto	1950-12-10	0.491	23.094	0.035
1960	1078512	De-facto	1960-12-11	0.493	22.818	0.030
1970	1428082	De-facto	1970-05-10	0.493	20.927	0.023
1980	1795012	De-facto	1980-05-11	0.495	24.146	0.026
1990	2330789	De-facto	1990-05-13	0.494	25.890	0.020
2000	2839177	De-facto	2000-05-14	0.495	27.715	0.018
2010	3405813	De-facto	2010-05-16	0.497	30.014	0.014
2023	4064780	De-jure	2023-02-01	0.504	32.864	NA

Cabe aclarar que no nos adentraremos en detalle sobre la evolución y situación actual demográfica de Panamá, ni daremos estimaciones más acabadas que las publicadas por el INEC (por ejemplo solo tomamos población censada, no ajustada), sino que nos servirá para trabajar en R distintas aplicaciones.

Edad y sexo

Visualicemos su distribución por edad y sexo. Lo primero que se observa es un proceso de envejecimiento, con angostamiento en la base, y resalta a simple vista una declaración por edad más precisa con el tiempo que luego mediremos más en detalle.

panama %>%
     filter(!Año %in% c(1970,1980), Edad < 100) %>% # el censo de 1970 se encuentra agrupado por edad quinquenal
     mutate(Pop = ifelse(Sexo == "Hombre", -Pop, Pop)/sum(Pop), .by = c(Año)) %>% 
     ggplot(aes(Edad, Pop, fill = Sexo)) +
     geom_line() +
     scale_x_continuous(labels = seq(0,100,10), breaks =  seq(0,100,10)) +
     scale_y_continuous(labels = abs) +
     coord_flip() +
     ggtitle("Pirámides de población censal. Panamá") +
     theme_bw() + 
     facet_wrap(~Año)

El último censo nos muestra una base más pequeña, sugiriendo una evolución intercensal de nacimientos decreciente, una posible omisión en menores, o ambos factores.

panama %>%
     filter(Año == 2023) %>%
     mutate(Pop = ifelse(Sexo == "Hombre", -Pop, Pop)/sum(Pop), .by = c(Año)) %>% 
     ggplot(aes(Edad, Pop, fill = Sexo)) +
     geom_bar(stat = "identity") + 
     scale_fill_manual(values = c("blue", "orange")) +
     scale_x_continuous(labels = seq(0,100,10), breaks =  seq(0,100,10), limits = c(0, 100)) +
     scale_y_continuous(labels = abs) +
     coord_flip() +
     theme_bw() + 
     ggtitle("Pirámides de población censal 2023. Panamá")

El índice de masculinidad permite ver el proceso de feminización visto en la tabla resumen, sobre todo en edades activas. El IM en censos previos a 1990 muestra una gran variabilidad, por lo que utilizamos la función geom_smooth de ggplot para tener idea de su patrón. Aquí puede tener una gran influencia el cambio en el tipo de censo, donde el último realizado es el primero de-jure (o de derecho), enfocándose en captar la población residente.

panama %>% 
     filter(Edad < 100, Año != 1970) %>% 
     summarise(IM = Pop[Sexo=="Hombre"]/Pop[Sexo=="Mujer"], .by = c(Edad, Año)) %>% 
     ggplot(aes(Edad, IM, col = factor(Año))) +
     geom_line(data = . %>% filter(Año < 2023)) +
     geom_smooth(data = . %>% filter(Año < 2023), span = .6, se = F) +
     geom_line(data = . %>% filter(Año == 2023), col = 1, linewidth = 1.5) +
     geom_hline(yintercept = 1, col = 2 , linetype = 2) +
     labs(color = "Censo") +
     scale_x_continuous(labels = seq(0, 100, 5), breaks = seq(0, 100, 5)) +
     ggtitle("Índice de masculinidad (Hombre/Mujer) por edad. Censos de Panamá") +
     theme_bw()

En términos de análisis intercensal por cohorte, en ausencia de importantes eventos de impacto demográfico (flujos migratorios, guerra, desplazamientos, etc.), y en ausencia también de omisión diferencial por censo, se esperaría ver que el ratio por cohorte refleje el efecto de la mortalidad (Siegel and Swanson 2004). No puede dar idea de omisión diferencial sino de consistencia. Su cálculo:

\[ cr_x(t, t+n) = \frac{P_{x+n}(t+n)}{P_{x}(t)} \]

El ratio se puede calcular por edad simple o por grupos quinquenales para evitar fluctuaciones. En los dos últimos períodos intercensales tenemos distintos intervalos, pero igual nos permite conocer algo. Asumamos que la distancia es de 13 años exactos en el último período intercensal:

panama_2010_2023_cohort_ratios <- panama %>% 
     filter(Año %in% c(2010, 2023), Edad < 80) %>% 
     # distinguimos cada cohorte
     mutate(cohorte = Año - Edad - 1) %>%
     # Sumamos población por cohorte y la dividimos
     summarise(cr = Pop[Año == 2023]/Pop[Año == 2010], 
               .by = c(cohorte, Sexo)) %>% 
     arrange(-cohorte)
head(panama_2010_2023_cohort_ratios)

##   cohorte   Sexo        cr
## 1    2009 Hombre 1.0122336
## 2    2009  Mujer 0.9957262
## 3    2008 Hombre 1.1131003
## 4    2008  Mujer 1.1198611
## 5    2007 Hombre 1.0344725
## 6    2007  Mujer 1.0188749

Aquellos valores por encima de 1 pueden hablar de omisión diferencial entre censos o fenómenos demográficos específicos (aunque de vuelta, el método ve consistencia). Por ejemplo la población nacida entre 2005 y 2009 parece estar subenumerada en 2010 (la población censada 13 años después es mayor). Luego en las cohortes 1995-1980 (con edades activas en período intercensal) se observa un ratio mayor a 1 en mujeres, y en varones cierta estabilidad. Luego en cohortes más antiguas el comportamiento es decreciente, como se espera. Ojo! esta condicionado al cambio de tipo de censo también.

panama_2010_2023_cohort_ratios %>%
     filter(cohorte > 1950) %>% 
     ggplot(aes(cohorte, cr, col = Sexo)) +
     geom_hline(yintercept = 1, linetype = 2) +
     geom_line() +
     scale_x_continuous(trans = "reverse", labels = seq(1920,2010,10), breaks = seq(1920,2010,10)) +
     scale_y_continuous(labels = seq(0,1.5,.1), breaks = seq(0,1.5,.1)) +
     scale_color_manual(values = c("blue", "orange")) +
     ggtitle("Relación de cohortes por edad simple, según sexo. Panamá 2023") +
     theme_bw()

En algunos casos, una visión por grupos quinquenales de edad puede ser más ilustrativa visualmente, removiendo la aleatoriedad propia de la desagregación por edad simple. En base a estos primeros diagnósticos visuales podemos identificar algunos aspectos que podemos analizar o ajustar mediante técnicas demográficas específicas.

Edad ignorada

Empecemos por algo bastante común: los casos de edad y/o sexo ignorada. El año 2023 tiene casos de edad ignorada pero no sexo ignorado.

panama_2023 <- panama %>% 
     filter(Año == 2023) %>% 
     mutate(Edad = pmin(Edad, 100)) %>% 
     summarise(Pop = sum(Pop), .by = c(Sexo, Edad))

# filtro casos con datos vacíos
panama_2023 %>% filter(is.na(Edad) | is.na(Sexo))

##     Sexo Edad Pop
## 1 Hombre   NA 291
## 2  Mujer   NA  44

Veamos el impacto:

# cuantas personas de edad ignorada
sum(panama_2023$Pop[is.na(panama_2023$Edad)])

## [1] 335

# porcentaje en el total censado
sum(panama_2023$Pop[is.na(panama_2023$Edad)])/sum(panama_2023$Pop)*100

## [1] 0.008241528

El impacto en este caso es menor, por lo que la redistribución proporcional es un buen recurso. La idea es tomar los casos de edad desconocida, obtener su suma por sexo, y mediante la proporción con edad y sexo dentro de lo conocido, imputar. Luego debemos checkear que tenemos el mismo total poblacional.

# redistribución
panama_2023_red <- panama_2023 %>% 
     filter(is.na(Edad)) %>% 
     summarise(Pop_unk = sum(Pop), .by = Sexo) %>% 
     right_join(panama_2023 %>% 
                     filter(!is.na(Edad))) %>% 
     mutate(Pop_red = Pop + Pop/sum(Pop) * Pop_unk, .by = Sexo)

# chequeemos que el total se respete
sum(panama_2023$Pop) == sum(panama_2023_red$Pop_red)

## [1] TRUE

Renombremos para continuar análisis:

panama_2023 <- panama_2023_red %>% mutate(Pop = Pop_red) %>% select(-Pop_red, -Pop_unk)

En el caso de que se tenga más de una variable con desconocido, típicamente sexo, el paso previo debe ser repetido: una opción es primero distribuir el ignorado de edad dentro de cada sexo, y luego el sexo.

Atracción de edad y suavizamiento

La distribución por edad puede mostrar algunos picos que podrían atribuirse a sesgos en la declaración, por atracción hacia ciertos dígitos, o hacia edades más avanzadas. Existen varios indicadores de atracción implementados en DemoTools para su diagnóstico. Por ejemplo para el año 2023, en hombres, solo para edades entre 10 y 65:

library(DemoTools)
panama_2023_hombre <- panama_2023 %>% filter(Sexo == "Hombre") 

# Whipple
whipple_panama2023_hombre <- check_heaping_whipple(panama_2023_hombre$Pop, panama_2023_hombre$Edad, 
                             ageMin = 10, ageMax = 65, digit = c(0, 5))
# Myers
myers_panama2023_hombre <- check_heaping_myers(panama_2023_hombre$Pop, panama_2023_hombre$Edad, 
                             ageMin = 10, ageMax = 65)
# Bachi
bachi_panama2023_hombre <- check_heaping_bachi(panama_2023_hombre$Pop, panama_2023_hombre$Edad, 
                             ageMin = 10, ageMax = 65)
# juntar
data.frame(Indicador = c("Whipple", "Myers", "Bachi"),
           Valor = c(whipple_panama2023_hombre, myers_panama2023_hombre, bachi_panama2023_hombre)) %>% 
     kbl(caption = "Indicadores de atracción de edades. Panamá 2023") %>% 
     kable_classic(full_width = F)

Indicadores de atracción de edades. Panamá 2023
Indicador	Valor
Whipple	1.0122304
Myers	0.6557403
Bachi	0.5601954

Respecto al rango de valores “aceptables”, Whipple (con resultados para dígitos 0 y 5) va de 1 (sin heaping) a 5 (concentración total), Myers mide la atracción para todos los dígitos y va de 0 (sin heaping) a 90 (mayor atracción), siendo Bachi muy similar a este último. Como se ve, la declaración se puede considerar muy buena. Estos indicadores de preferencia/atracción descansan en una comparación respecto a un supuesto de distribución por edad (en general lineal), donde a veces la estructura por edad puede ser resultado de un efecto cohorte.

Quizás convenga también verlo desde una perspectiva histórica: el reciente censo presenta una declaración mucho más precisa que los anteriores. También se ve una alta correlación entre índices.

panama_heaping_index <- panama %>% 
     # dejamos fuera Edad ignorada y censo 1970 por edades agrupadas
     filter(!is.na(Edad), Año != 1970) %>% 
     # dividir en listas y aplicar funciones de cálculo de indicadores
     split(list(.$Año, .$Sexo)) %>% 
     map_df(function(census){
          whipple <- check_heaping_whipple(census$Pop, census$Edad, ageMin = 10, ageMax = 65)
          myers <- check_heaping_myers(census$Pop, census$Edad, ageMin = 10, ageMax = 65)
          bachi <- check_heaping_bachi(census$Pop, census$Edad, ageMin = 10, ageMax = 65)
          data.frame(Indicador = c("Whipple", "Myers", "Bachi"),
                     Valor = c(whipple, myers, bachi),
                     Año = unique(census$Año),
                     Sexo = unique(census$Sexo))
     })

panama_heaping_index %>% 
     ggplot(aes(Año, Valor, col = Sexo))+
     geom_line() + geom_point() +
     facet_wrap(~Indicador, scales = "free") +
     ggtitle("Evolución censal de indicadores de atracción por edad según sexo. Panamá") +
     scale_color_manual(values = c("blue", "orange")) +
     theme_bw()

En el protocolo de WPP-NU (United Nations and Social Affairs 2022), la División de Población consideró el índice de Bachi para los casos de edad simple, junto con el nivel medio de educación completado por los adultos de la población, para determinar el grado de suavizamiento sobre la edad a aplicar mediante una media móvil.

La población por edad en el último censo tiene la siguiente distribución:

panama_2023 %>%
     ggplot(aes(Edad, Pop, col = Sexo)) +
     geom_point() +
     scale_x_continuous(labels = seq(0, 100, 5), breaks = seq(0, 100, 5)) +
     theme_bw()+
     scale_color_manual(values = c("blue", "orange")) +
     ggtitle("Población censal por edad y sexo. Panamá 2023") +
     theme(panel.grid.minor.x = element_blank())

Lo ideal sería contar con una encuesta post-censal que permita inferir la distribución por edad más acorde. Cuando no es posible, distintas opciones están disponibles para suavizar. Aquí vemos algunas implementadas en DemoTools, lascuales funcionan sobre edades agrupadas, por lo que ese paso es necesario previamente:

graduación <- panama_2023 %>% 
     # separamos los dos sexos en una lista y aplicamos la graduación según cada método
     split(.$Sexo) %>% 
     map_df(function(panama_2023_sexo){
          # agrupamiento en edades quinquenales
          panama_2023_sexo_5 <- panama_2023_sexo %>% 
               mutate(Edad5 = trunc(Edad/5) * 5) %>% 
               summarise(Pop5 = sum(Pop), .by = Edad5)
          # cómputo de indicadores
          sprague <-   graduate(panama_2023_sexo_5$Pop5, panama_2023_sexo_5$Edad5, method = "sprague")
          beers <-     graduate(panama_2023_sexo_5$Pop5, panama_2023_sexo_5$Edad5, method = "beers(ord)")
          mono.grad <- graduate(panama_2023_sexo_5$Pop5, panama_2023_sexo_5$Edad5, method = "mono")
          pclm.grad <- graduate(panama_2023_sexo_5$Pop5, panama_2023_sexo_5$Edad5, method = "pclm")
          # resulatdo por sexo
          results_sex <- data.frame(Sexo = unique(panama_2023_sexo$Sexo),
                                    Edad = rep(0:100, 4),
                                    Pop = c(sprague, beers, mono.grad, pclm.grad),
                                    method = c(rep("sprague", 101), rep("beers", 101), 
                                               rep("mono", 101), rep("pclm", 101)))
          return(results_sex)
     })

# plot
plot_graduación <- bind_rows(panama_2023 %>% select(Sexo, Edad, Pop) %>% mutate(method = "Observado"),
          graduación) %>% 
     ggplot() +
     geom_point(data = . %>% filter(method == "Observado"), aes(Edad, Pop), size = .5) +
     geom_line(data = . %>% filter(method != "Observado"), aes(Edad, Pop, col = method), linewidth = 1) +
     theme_bw() +
     ggtitle("Población censal por edad y sexo, y métodos de suavizamiento. Panamá 2023") +
     facet_wrap(~Sexo)
plot_graduación

Excepto en edades finales, la performance es bastante similar entre métodos. También se puede utilizar métodos más simples como promedios móviles, en este caso con una ventana de 3 años utilizando la función rollmean del paquete zoo para edades desde 5:

mav <-panama_2023 %>% 
     split(.$Sexo) %>% 
     map_df(function(panama_2023_sexo){
     mav <- c(panama_2023_sexo$Pop[panama_2023_sexo$Edad < 5],
                   zoo::rollmean(panama_2023_sexo$Pop[panama_2023_sexo$Edad >= 5], 
                                 k = 3, fill = "extend"))  
     results_sex <- data.frame(Sexo = unique(panama_2023_sexo$Sexo),
                               Edad = panama_2023_sexo$Edad,
                               Pop = mav,
                               method = "mav")
     return(results_sex)
     })
          
plot_graduación +
     geom_line(data = mav, aes(Edad, Pop, col = method), linewidth = 1)

Siempre es importante corroborar totales (y distribución por grandes grupos de edad), aunque algunos métodos lo garanticen de antemano (por ejemplo Sprague):

graduación %>% 
     summarise(Pop_grad = sum(Pop), .by = method) %>% 
     mutate(Diff = sum(panama_2023$Pop) - Pop_grad)

##    method Pop_grad        Diff
## 1 sprague  4064780  0.00000000
## 2   beers  4064780  0.00000000
## 3    mono  4064780  0.00000000
## 4    pclm  4064780 -0.01501837

En el caso de WPP-NU (United Nations and Social Affairs 2022), las evaluaciones de acumulación de edad y el suavizado se realizaron por separado para las edades de niños y adultos y las series resultantes se sumaron al postajuste. Cuando la serie de datos del censo se agrupó por grupos de edad de cinco años, la serie ajustada y suavizada finalmente se graduó a edades simples utilizando spline monótono (method = "mono"). Cuando fue necesario, los recuentos se redistribuyeron al grupo de edad indefinido de 105 años o más utilizando la función OPAG (documentación), que aplicó la tabla de vida estimada para el país y el año del censo como un estado estable.

Edades iniciales

Vimos en la pirámide que las edades iniciales sugieren una reciente baja de la natalidad, sin embargo el descenso de la natalidad parece comenzar en 2019. Evaluaremos la población en edades menores a 10 sobreviviendo los nacimientos de esas cohortes.

panama_nacs <- births_LA %>% 
     # distinct %>% 
     filter(País == "Panama", Edad_interv != 1) %>% 
     group_by(País, Año, Sexo, Edad) %>% slice(1) %>% ungroup %>%  
     summarise(Births = sum(Births), .by = c(Año, Sexo))
panama_nacs %>% 
     ggplot(aes(Año, Births, col = Sexo)) +
     geom_line() + geom_point() +
     ggtitle("Nacimientos por sexo según año de ocurrencia. Panamá") +
     scale_color_manual(values = c("blue", "orange")) +
     theme_bw()

Por el lado de las defunciones, nos interesan aquellas de edad menor a 10 desde el año 2013:

panama_defs <- deaths_LA %>% 
     # distinct %>% 
     filter(País == "Panama", Año >= 2013) %>% 
     group_by(País, Año, Sexo, Edad, Edad_interv) %>% slice(1) %>%   
     summarise(Deaths = sum(Deaths)) %>% 
     arrange(Sexo, Edad, Año) %>% ungroup

El año 2022 es el único con edades quinquenales (podemos mirarlo). Utilicemos pclm para desagrupar:

panama_defs_2022 <- panama_defs %>% 
     filter(Año == 2022, !is.na(Edad)) %>% 
     split(.$Sexo) %>% 
     map_df(function(panama_defs_2022_sexo){
          panama_defs_2022_sexo_grad <- graduate(panama_defs_2022_sexo$Deaths,
                                                 panama_defs_2022_sexo$Edad, method = "pclm")
          return(data.frame(Sexo = unique(panama_defs_2022_sexo$Sexo), 
                            Edad = 0:85, 
                            Deaths = panama_defs_2022_sexo_grad))
     })

# incorporar
panama_defs <- 
     bind_rows(panama_defs %>% filter(Año != 2022),
               panama_defs_2022 %>% mutate(Año = 2022))

El siguiente paso es asignar el año de nacimiento en cada año y edad, es decir distribuir por triángulos inferior y superior. Para esto es necesario el concepto de \(a_x\): la cantidad de años promedio vivida por las defunciones de edad \(x\), y desde allí entender que este tiempo promedio y el porcentaje de defunciones correspondiente al triángulo superior tienen una relación de equivalencia formal (Andreev and Kingkade 2015). En el siguiente plot podemos ver visualmente de qué estamos hablando utilizando el paquete LexisPlotR:

library(LexisPlotR)
triang_inf <- data.frame(group = c(1, 1, 1),
                       x = c("2021-01-01", "2022-01-01", "2022-01-01"),
                       y = c(0, 0, 1))
triang_sup <- data.frame(group = c(1, 1, 1),
                       x = c("2021-01-01", "2021-01-01", "2022-01-01"),
                       y = c(0, 1, 1))
lexis_grid(year_start = 2020, year_end = 2023, age_start = 0, age_end = 2) %>% 
     lexis_polygon(x = triang_inf$x, y = triang_inf$y, group = triang_inf$group, fill = "yellow") %>%
     lexis_polygon(x = triang_sup$x, y = triang_sup$y, group = triang_sup$group, fill = "blue") +
     ggtitle("Triángulo superior e inferior de las defunciones de edad 0\n en el año 2021. Diagrama teórico") +
     annotate("text", x = as.Date("2021-04-01"), y = .75, label = "a(0)", col = "white") +
     annotate("text", x = as.Date("2021-10-01"), y = .25, label = "1-a(0)", col = "black")

Para obtener \(a_0\) existen muchas opciones. DemoTools pone a disposición las opciones de Coale-Demeny, la utilizada por el grupo HMD a través de la función lt_rule_1a0, o Keyfitz-Flieger que toman como argumento la tasa de mortalidad infantil (definida como \(D_0(t)/B(t)\)) (Preston, Heuveline, and Guillot 2000; Database 2021; Wachter 2006). Tomemos el promedio del período como referencia:

a0 <- panama_nacs %>% 
     inner_join(panama_defs %>% filter(Edad == 0)) %>% 
     mutate(a0 = lt_rule_1a0(rule = "ak", q0 = Deaths/Births)) %>%
     summarise(mean(a0)) %>% pull()

Para las edades 1:4 podemos utilizar los factores de separación de Glover y luego .5 (Nations 1983).

ax <- data.frame(age = 0:9, 
                 ax = c(a0, 0.41, 0.47, 0.48, 0.48, rep(.5, 5)))

Para restar las defunciones a cada cohorte, podemos aplicar la función lexis_pop (tomada del script funs.R), que brinda como output la población a cada año y un gráfico con el diagrama de Lexis.

# prepara data para lexis
panama_defs_lexis <- panama_defs %>% 
     filter(Edad %in% 0:9, Año %in% 2013:2022) %>% 
     select(yd = Año, age = Edad, sex = Sexo, D = Deaths) 
     # distinct()
panama_nacs_lexis <- panama_nacs %>% 
     filter(Año %in% 2013:2022) %>% 
     select(yb = Año, sex = Sexo, B = Births) 
     # distinct()

# lexis
pop_0a9_2023_mujeres <- lexis_pop(deaths = panama_defs_lexis %>% filter(sex == "Mujer"), 
                                  births = panama_nacs_lexis %>% filter(sex == "Mujer"), 
                                  t0 = 2013, ax = ax)
pop_0a9_2023_hombres <- lexis_pop(deaths = panama_defs_lexis %>% filter(sex == "Hombre"), 
                                  births = panama_nacs_lexis %>% filter(sex == "Hombre"), 
                                  t0 = 2013, ax = ax)

# por ejemplo para Hombres
pop_0a9_2023_mujeres$plot_lexis +
     ggtitle("Evolución por cohorte según eventos observados. \nDiagrama de Lexis. Panamá")

Finalmente, volvemos al objetivo final de esta sección, comparar la población menor a 10 con la censada en 2023. Tengamos en cuenta la diferencia del momento de evaluación (fecha de referencia censal y comienzo de año), aunque es solo de un mes. El gráfico siguiente y la comparación de las diferencias absoluta y relativas nos permite tener una idea más acabada sobre si un ajuste es necesario:

# población estimada por Lexis
panama_lexis_pop2023 <- bind_rows(
     pop_0a9_2023_mujeres$pop_out %>% mutate(Sexo = "Mujer"),
     pop_0a9_2023_hombres$pop_out %>% mutate(Sexo = "Hombre")) %>% 
     filter(t == 2023)

# censo
panama_2023 <- census_LA %>% 
     filter(País == "Panama", Año == 2023)

ggplot() +
     geom_line(data = panama_2023 %>% filter(Edad < 15), aes(Edad, Pop))+
     geom_line(data = panama_lexis_pop2023, aes(age, pop), col = 2) +
     scale_y_continuous(limits = c(0, 40000)) +
     scale_x_continuous(labels = 0:15, breaks = 0:15) +
     facet_wrap(~ Sexo) +
     ggtitle("Población censal y estimación a aprtir de eventos. Panamá 2023")  +
     theme_bw() +
     theme(panel.grid.minor.x = element_blank())

En términos cuantitativos, la diferencia es de alrededor del 7% en 0-4,mientras de 3% en 5-9.

# censo y estimación por lexis en edades 0-4 y 5-9
censo <- panama_2023 %>% 
     filter(Edad < 10) %>% 
     mutate(Edad = ifelse(Edad < 5, "0-4", "5-9")) %>% 
     summarise(Pop_censal = sum(Pop), .by = Edad) 
lexis <- panama_lexis_pop2023 %>% 
     filter(age < 10) %>% 
     mutate(Edad = ifelse(age < 5, "0-4", "5-9")) %>% 
     summarise(Pop_estimada = trunc(sum(pop)), .by = Edad) 

# tabla resumen
inner_join(censo, lexis) %>% 
     mutate(Dif_Abs = Pop_censal - Pop_estimada, 
            Dif_Rel = round((Pop_censal - Pop_estimada)/Pop_censal,2)) %>% 
     kbl(caption = "Estimación poblacional y censal de menores de 10 años, según edad. Panamá 2023 (inicio de año y fecha de referencia)") %>% 
     kable_classic(full_width = F)

Estimación poblacional y censal de menores de 10 años, según edad. Panamá 2023 (inicio de año y fecha de referencia)
Edad	Pop_censal	Pop_estimada	Dif_Abs	Dif_Rel
0-4	322449	344226	-21777	-0.07
5-9	357880	369801	-11921	-0.03

Mortalidad

Construcción y graduación

Obtengamos la tabla de mortalidad de Panamá para el año 2022. Para eso, requerimos una estimación de la población expuesta de mujeres por edades quinquenales (dado que las defunciones lo están), por lo que haremos una interpolación por grupo de edad y sexo para correr desde la fecha de referencia censal (“01/02/2023”) hacia mitad de año de 2022, según crecimiento 2010-2023. Atención! Utilizamos este supuesto de interpolación por simplificación, pero puede ser especialmente problemático en casos de alto impacto por covid19 si vamos hacia los años 2020 y 2021 (en general ante cualquier efecto fuerte de período). En Panamá se observó un 25% más de defunciones en 2020 respecto a años previos, por lo que la exposición por grupos de edad (denominador de tasas) puede variar significativamente por su concentración en edades adultas y finales y por su estacionalidad.

# población al 2010 y 2023
panama_2010_2023 <- panama %>%
     filter(Año %in% c(2010, 2023)) %>% 
     mutate(Date = dec.date(dmy(Fecha)),
            Edad = pmin(Edad, 85, na.rm = T)) %>% 
     summarise(Pop = sum(Pop), .by = c(Date, Sexo, Edad))

# estimar exposición al 2022
panama_exposures_2022 <- bind_rows(
     interp_tidy(panama_2010_2023 %>% filter(Sexo == "Mujer"), Edad, Date, Pop, 2022.5) %>% 
          mutate(Sexo = "Mujer"),
     interp_tidy(panama_2010_2023 %>% filter(Sexo == "Hombre"), Edad, Date, Pop, 2022.5) %>% 
          mutate(Sexo = "Hombre"))

# agrupación en edades quinquenales
panama_exposures5_2022 <- panama_exposures_2022 %>% 
     mutate(Edad = ifelse(Edad %in% 1:4, 1, pmin(trunc(Edad/5)*5, 85))) %>% 
     summarise(Pop = sum(Pop), .by = c(Sexo, Edad))

# obtener defunciones del año 2022
panama_deaths_mujer_2022 <- deaths_LA %>% 
     filter(País ==  "Panama", Año == 2022, Sexo == "Mujer") %>% 
     select(Edad, Deaths)

Finalmente obtenemos las tasas por edad

# tasas por edad para el caso de Mujer
panama_mx5 <- panama_exposures5_2022 %>% 
     filter(Sexo == "Mujer") %>% 
     inner_join(panama_deaths_mujer_2022) %>% 
     mutate(mx = Deaths/Pop,
            Edad_media = Edad + 2.5)

panama_mx5_plot <- panama_mx5 %>% 
     ggplot(aes(Edad, mx)) + 
     geom_step() +
     scale_y_log10()+
      ggtitle("Tasa de mortalidad por edades agrupadas a partir de \ndatos censales y defunciones observadas. Panamá 2022") +
     theme_bw() 
panama_mx5_plot

Si bien no lo incorporamos, el hecho de incluir o no el ajuste a la población menor a 10 impacta en las tasas, disminuyéndolas. Si nos interesa extender o ajustar las tasas en edades finales podemos hacer uso del paquete MortalitiLaws, para ajustar y predecir tasas en edades finales. Probemos con las leyes makeham y kannisto, estimando sus parámetros en base a las tasas de las edades 50 a 80, y extendiendo hacia 100 la edad del grupo abierto final.

library(MortalityLaws)

# distintas leyes de ajuste
makeham_model <- MortalityLaw(x = panama_mx5$Edad, mx = panama_mx5$mx, law = "makeham", fit.this.x = seq(50, 80, 5)) %>% predict(x = seq(50, 100, 5))
kannisto_model <- MortalityLaw(x = panama_mx5$Edad, mx = panama_mx5$mx, law = "kannisto", fit.this.x = seq(50, 80, 5)) %>% predict(x = seq(50, 100, 5))
panama_mx5_laws <- data.frame(Edad = rep(seq(50, 100, 5), 2) + 2.5,
                              nMx = c(makeham_model, kannisto_model),
                              Model = c(rep("Makeham", 11), rep("Kannisto", 11)))
panama_mx5_plot +
     geom_point(data = panama_mx5_laws,
                aes(Edad, nMx, col = Model)) +
     ggtitle("Tasa de mortalidad por edades agrupadas a partir de \ndatos censales y defunciones observadas, y extensión hacia edades finales. Panamá 2022")

Como se ve, ambas estimaciones se diferencian a partir de la edad 80, ya que desde la edad 50 se utilizó para su ajuste. Por suerte,esta funcionalidad esta ya presente en las funciones de DemoTools, para crear una tabla de mortalidad. Seleccionemos Makeham.

# construcción
lt_panama_2022 <- lt_abridged(nMx = panama_mx5$mx, Age = panama_mx5$Edad, Sex = "f", OAnew = 100, extrapLaw = "makeham", extrapFit = seq(50, 80,5))
head(lt_panama_2022)

##    Age AgeInt          nMx       nAx         nqx        lx       ndx       nLx
## 0    0      1 0.0128969178 0.1225234 0.012752600 100000.00 1275.2600  98880.99
## 1    1      4 0.0008132085 1.5032515 0.003246243  98724.74  320.4845 394098.79
## 5    5      5 0.0002572245 2.5000000 0.001285296  98404.26  126.4786 491705.08
## 10  10      5 0.0003435139 2.5000000 0.001716096  98277.78  168.6541 490967.25
## 15  15      5 0.0004447108 2.5000000 0.002221085  98109.12  217.9087 490000.84
## 20  20      5 0.0006137075 2.5000000 0.003063837  97891.21  299.9227 488706.26
##           Sx      Tx       ex
## 0  0.9859596 8152030 81.52030
## 1  0.9974143 8053149 81.57174
## 5  0.9984994 7659051 77.83251
## 10 0.9980316 7167345 72.92946
## 15 0.9973580 6676378 68.05053
## 20 0.9967735 6186377 63.19645

También podemos graduar hacia edades simples ya que la función lt_abridged2single:

# graduación a edades simples
lt_panama_2022_simple <- lt_abridged2single(nMx = panama_mx5$mx, Age = panama_mx5$Edad, Sex = "f", OAnew = 100, extrapLaw = "makeham", extrapFit = seq(50, 80, 5))
panama_mx5_plot +
     geom_line(data =lt_panama_2022_simple, aes(Age, nMx)) +
     ggtitle("Tasa de mortalidad por edades agrupadas a partir de \ndatos censales y defunciones observadas, y graduación por edad simple. Panamá 2022")

Comparemos la \(e_0\) entre el cálculo por edades agrupadas y simples, y con otras fuentes, como por ejemplo WPP.

lt_panama_2022$ex[lt_panama_2022_simple$Age==0]

## [1] 81.5203

lt_panama_2022_simple$ex[lt_panama_2022_simple$Age==0]

## [1] 81.00496

Uso de tablas modelo

Existen situaciones donde nuestra confianza en las estadísticas vitales no es tan buena, y solo contamos con indicadores resumen a partir de estimaciones indirectas como \(q_0\), \({}_{5}q_0\), \(e_0\) o \({}_{45}q_{15}\). O también nos puede interesar saber qué patrón por edad seguiría la mortalidad de nuestra población bajo una de las familias de tablas de mortalidad (UN?). Para esto tenemos dos funciones implementadas por Sara Hertog, adaptando funciones del en software MortPack (disponible online aquí para aquellos más amigos de excel).

Por ejemplo bajo la familia Oeste de Coale Demeny, teniendo como insumo solo \(e_0\), esta función encuentra, dado el nivel, la interpolación de tasas que cumple con el valor dado:

e0 <- lt_panama_2022_simple$ex[lt_panama_2022_simple$Age==0]
panama_2022_mlt_one <- lt_model_cdun_match_single(type = "CD_West",
                           indicator = "e0", 
                           value = e0,
                           Sex = "f", OAnew = 100)

Pero también podemos tener dos valores, y buscar el nivel que mejor ajuste dentro de una familia dada. Intentémos con \({}_1q_0\) y \({}_45q_{15}\):

q0_1 <- lt_panama_2022_simple$nqx[lt_panama_2022_simple$Age==0]
q15_45 <- 1-lt_panama_2022_simple$lx[lt_panama_2022_simple$Age==60]/lt_panama_2022_simple$lx[lt_panama_2022_simple$Age==15]
panama_2022_mlt_two <- lt_model_cdun_combin_single(type = "CD_West", 
                            Sex = "f", 
                            q1 = q0_1, 
                            q5 = NA, 
                            indicator_adult = "45q15", 
                            value_adult = q15_45, OAnew = 100)

Veamos visualmente el resultado de ambos intentos:

panama_2022_mlts <- bind_rows(lt_panama_2022_simple %>% mutate(MLT = "Observed"),
                              panama_2022_mlt_one %>% mutate(MLT = "CD_West_e0"),
                              panama_2022_mlt_two %>% mutate(MLT = "CD_West_q0_&_45q15"))
panama_2022_mlts %>% 
     ggplot() +
     geom_point(data = . %>% filter(MLT == "Observed"), aes(Age, nMx)) + 
     geom_point(data = . %>% filter(MLT != "Observed"), aes(Age, nMx, col = MLT)) +
     scale_y_log10()+
     scale_x_continuous(breaks = seq(0,100,10), labels =  seq(0,100,10), limits = c(0,100)) +
     ggtitle("Tasa de mortalidad por edades simples a partir de \ndatos censales y defunciones observadas, \ny tablas modelo CD familia West. Panamá 2022") +
     theme_bw()

¿Qué conclusiones podemos obtener?

Otra opción interesante es la de utilizar el modelo log-quad (Wilmoth et al. 2012), para saber, bajo el patrón de experiencia histórica HMD, cuál serían las tasas por edad que corresponderían a nuestro valor de \({}_1q_0\). Esta función esta implementada en paquetes como MortCast o MortalityEstimate, pero también en DemoTools!

Fecundidad

Graduación

Calculemos las tasas abreviadas de fecundidad por edad para Panamá en el año 2022. Para eso, utilicemos la estimación de la población expuesta de mujeres hecha previamente, y recuperemos los nacimientos por edad de la madre.

# obtener nacimientos en 2022
panama_births5_2022 <- births_LA %>% 
     filter(País ==  "Panama", Año == 2022, Edad %in% 15:45) %>% 
     summarise(Births = sum(Births), .by = c(Edad))

# calcular tasas
panama_asfr5 <- panama_exposures5_2022 %>% 
     filter(Sexo == "Mujer") %>% 
     inner_join(panama_births5_2022) %>% 
     mutate(asfr = Births/Pop,
            Edad_media = Edad + 2.5)
panama_asfr5 %>% 
     ggplot(aes(Edad, asfr)) +
     geom_step() +
     ggtitle("Tasas de fecundidad por edades quinquenales a partir de \ndatos censales y nacimientos observados. Panama 2022") +
     theme_bw()

Aquí trataremos los métodos de splines calibrados de (Schmertmann 2014) (código disponible aquí), Beers (misma método utilizado en secciones previas), y el método de optimización cuadrática (quadratic optimization) (Grigoriev et al. 2018). Diversas pruebas hechas por el grupo de HFD muestran que en términos de flexibilidad y ajuste, el orden decreciente sugerido es CS y Beers. Intentemos qué resulta para nuestro caso:

# Beers
asfr_beers <-  graduate_beers(Value = panama_asfr5$asfr, Age = panama_asfr5$Edad, 
                              method = 'mod', OAG = FALSE ) * 5
# CS
asfr_cs <- graduate_cs(Value = panama_asfr5$asfr, Age = panama_asfr5$Edad)

# quadratic optimization
library("quadprog")
asfr_qo <- QOSplit(c(panama_asfr5$asfr), L=seq(15,45,5), AgeInt=rep(5,7))

# juntar ambas
asfr_graduation <- data.frame(
     Edad = rep(15:49, 3),
     Method = c(rep("Beers", 35), rep("CS", 35), rep("QO", 35)),
     asfr = c(asfr_beers, asfr_cs, asfr_qo$ASFR))
asfr_graduation %>% 
     ggplot(aes(Edad, asfr)) +
     geom_line(aes(col = Method)) +
     geom_point(data = panama_asfr5, aes(Edad_media, asfr)) +
     ggtitle("Tasas de fecundidad por edades quinquenales y \ngraduación por edades simples. Panamá 2022") +
     theme_bw()

Finalmente, corroboremos la TGF y edad media resultantes:

asfr_graduation %>% 
     summarise(TFR = sum(asfr), 
               Edad_Media = sum(asfr*(Edad+.5))/sum(asfr),
               .by = Method) %>% 
     bind_rows(data.frame(Method = "Observed", TFR = sum(panama_asfr5$asfr) * 5, Edad_Media = sum(panama_asfr5$asfr*panama_asfr5$Edad_media) / sum(panama_asfr5$asfr)))

##     Method      TFR Edad_Media
## 1    Beers 2.065389   27.52279
## 2       CS 2.064146   27.61072
## 3       QO 2.065389   27.50685
## 4 Observed 2.065389   27.48852

Se puede incorporar 10-14 al ajuste. En algunos casos, para CS es mejor re-calibrar los coeficientes de splines para tener un buen ajuste según el aptrón subyacente que se utiliza (United Nations and Social Affairs 2022).

Otros

Procedimiento de ajuste proporcional iterativo

En diversas situaciones tenemos una distribución de stock/eventos demográficos según categorías que por provenir de fuentes distintas, sus marginales no respetan totales deseados. Casos típicos son la consistencia flujos migratorios bilaterales (origen-destino), o matrices de cambio de estado entre dos momentos del tiempo, la consistencia por edad de proyecciones en áreas menores respecto a áreas mayores, o absorber categorías desconocidas, entre muchas otras. La técnica que veremos aquí se la conoce como “Procedimiento de ajuste proporcional iterativo” (“Iterative Proportional Fitting”), ”Raking“, o “tabla cuadrada” en ámbitos más demogrpaficos. Es una herramienta simple muy utilizada por muchos Institutos de Estadística Nacionales. Es válido para celdas positivas (0 puede salvarse con 0,0001) y marginales no nulos. Existe una generalización a \(n\) variables mediante un enfoque estadístico de regresión log-linear. Como desventaja, no tiene una interpretación demográfica clara.

En este ejemplo veremos un caso donde debemos reportar la población por edad proyectada a 2025 para las áreas menores (Comunas) de la Ciudad Autónoma de Buenos Aires (CABA), y como insumos se tiene:

la población por edad proyectada de CABA a 2025,
la población total proyectada por Comuna a 2025,
la población censal por edad según Comuna a 2010.

Tenemos una distribución inicial (c), que en su suma de columnas (marginal) debe respetar el total por comuna (b), y la suma por fila debe respetar el total del área mayor por edad (a). Veamos la tabla:

proy_area_menor_totales <- read_xlsx("argentina_caba.xlsx", sheet = "proy_area_menor")
proy_area_mayor_edades <- read_xlsx("argentina_caba.xlsx", sheet = "proy_area_mayor")
censo_area_menor_edades_2010 <- read_xlsx("argentina_caba.xlsx", sheet = "censo_area_menor")
a <- proy_area_mayor_edades %>% filter(Anio == 2025) %>% 
     mutate(Edad = ifelse(Edad > 95, 95, Edad)) %>% 
     summarise(Pop = sum(Pop), .by = Edad) %>% 
     pull(Pop)
b <- proy_area_menor_totales %>% filter(Anio == 2025) %>% pull(Pop)
c <- censo_area_menor_edades_2010 %>% 
     select(-Edad_interv) %>% 
     pivot_wider(names_from = Comuna, values_from = Pop) %>%
     select(-Edad) %>% 
     as.matrix()

Si quisiéramos visualizar la matriz, tendríamos lo siguiente, donde los marginales están coloreados:

M <- data.frame(rbind(cbind(c, a), c(b, sum(b))))
colnames(M) <- c(1:15, "a")
rownames(M) <- c(seq(0, 95,5), "b")
kable(M, caption = "Matriz con datos para ajuste de marginales: población censal y proyecciones. Argentina") %>% 
     kable_classic(full_width = F) %>% 
     row_spec(21, background = "yellow") %>% 
     column_spec(column = 17, background = "lightblue")

Matriz con datos para ajuste de marginales: población censal y proyecciones. Argentina
	1	2	3	4	5	6	7	8	9	10	11	12	13	14	15	a
0	11669	6220	9889	15589	9006	9095	14368	16754	10091	9242	9899	11321	12023	10662	9810	194379
5	10149	5696	9227	15052	8659	8409	13633	15552	10195	9403	10424	10476	10650	9254	9593	196125
10	9913	5486	8986	15106	8295	7632	13200	14927	9985	8986	10043	9732	10010	8851	9349	195294
15	11582	8239	10211	15828	9390	8565	14038	15445	10444	9813	10730	10600	11288	11247	10261	196262
20	19494	17054	16523	18189	13298	11646	17561	17423	11672	11668	12499	13028	15555	19255	13260	196502
25	21421	16686	17794	18214	15948	14191	17836	16082	11934	12239	13971	15643	18736	21677	15222	198334
30	18893	13106	16784	17928	16469	15769	17755	14681	12341	12943	15131	18288	20630	20768	16583	209424
35	15667	9916	14233	15750	13902	13967	15749	13003	11699	11947	14080	15960	18468	16958	14027	223931
40	12899	8268	11898	13732	11416	11248	13484	11124	10341	10367	12380	12777	14944	14185	11813	220349
45	12043	8362	11179	12349	10752	10869	12772	9747	9652	10472	12095	12480	14078	13182	11594	216824
50	11991	8897	10769	11955	10948	11379	12858	8682	9493	10705	12305	12566	14122	13252	11099	197400
55	10972	9246	10431	10918	10372	10767	11800	7567	8913	9967	11414	11473	13613	13219	10464	162539
60	10069	9491	9585	9515	9750	10271	10933	6829	8341	9177	10711	11111	13651	12902	9779	149043
65	8327	8394	8095	7950	8231	8521	9114	5716	7156	7513	9019	9368	11799	11257	7955	142526
70	6437	6823	6620	6531	6723	7007	7582	4795	6003	6511	7391	7686	9402	9060	6602	125159
75	5497	5701	5934	5686	5937	6210	6883	4081	5567	5945	6817	6945	8396	7739	5958	104476
80	4512	4965	5071	4479	5209	5534	5864	2849	4420	4964	5796	5677	7197	6569	4964	75803
85	2811	3350	2876	2426	3132	3305	3451	1415	2522	2903	3562	3415	4513	3970	2948	46228
90	1173	1511	1072	835	1218	1296	1331	457	811	1013	1232	1208	1697	1493	1001	23434
95	367	521	360	213	350	395	379	108	217	244	333	362	559	470	292	12648
b	259205	149389	193884	241193	188053	186028	242921	230745	171921	170896	190225	215505	236781	227360	182574	3086680

Utilicemos el paquete mipfp, especializado en esto. El argumento seed es la matriz inicial, que debe ser ajustada, mientras que en target.list se referencia como 1 al marginal de fila, y como 2 al de columna.

# install.packages("mipfp")
library(mipfp)
censo_area_menor_edades_2025 <- Ipfp(seed = c, 
                                     target.list = list(1, 2), 
                                     target.data = list(a, b))

# check
round(colSums(censo_area_menor_edades_2025$x.hat) - b, 5)

##  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 
##  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0

round(colSums(t(censo_area_menor_edades_2025$x.hat)) - a, 5)

##  [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

El paquete permite incrementar las dimensiones de análisis si se desea. Un excelente tutorial referido a migraciones puede encontrarse en este Training Session: Tools for working with migration data in R (IPC2021).

Patrones de migración neta por edad

Veamos las funciones mig_un_fam (patrón dado según familias UN) y mig_estimate_rc (estimación bayesiana de parámetros de modelo Roger-Castro, tutorial aquí).

Actividad

Selecciona un país que te interese y obtén el gráfico de pirámides censales. ¿Qué llama más la atención?
Obtén la evolución temporal de 1 o más indicadores de atracción de dígitos en la edad.
Para algún país con datos de defunciones en algún año censal disponible, obtener las tasas de mortalidad por edad y construir una tabla de mortalidad mediante DemoTools. Probar distintos ajustes en edades finales.

Por último, si sos tan amable, puedes contestar esta breve encuesta? Gracias!

Referencias

Andreev, Evgeny M., and W. Ward Kingkade. 2015. “Average age at death in infancy and infant mortality level: Reconsidering the Coale-Demeny formulas at current levels of low mortality.” Demographic Research 33 (August): 363–90. https://www.demographic-research.org/articles/volume/33/13.

Database, Human Mortality. 2021. “Methods Protocol for the Human Mortality Database, Version 6.” https://www.mortality.org/File/GetDocument/Public/Docs/MethodsProtocolV6.pdf.

Estadística y Censos, Insituto Nacional de. 2016. “El Proceso de Transición Demográfica en Panamá” Sección 221. https://www.inec.gob.pa/archivos/P7441El%20Proceso%20de%20Transici%C3%B3n%20Demogr%C3%A1fica%20en%20Panam%C3%A1.pdf.

Grigoriev, Pavel, Anatoli I. Michalski, Vasily P. Gorlishchev, Dmitry Jdanov, and Vladimir M. Shkolnikov. 2018. “New methods for estimating detailed fertility schedules from abridged data.” ResearchGate, January. https://www.researchgate.net/publication/322266948_New_methods_for_estimating_detailed_fertility_schedules_from_abridged_data.

Nations, United. 1983. “Manual X: Indirect techniques for demographic estimation.” https://unstats.un.org/unsd/demographic/standmeth/handbooks/Manual%20X.pdf.

Preston, Samuel, Patrick Heuveline, and Michel Guillot. 2000. Demography: Measuring and Modeling Population Processes. Hoboken, NJ, USA: Wiley. https://www.wiley.com/en-us/Demography%3A+Measuring+and+Modeling+Population+Processes-p-9781557864512.

Schmertmann, Carl P. 2014. “Calibrated spline estimation of detailed fertility schedules from abridged data\({^1}\).” Rev. Bras. Estud. Popul. 31 (December): 291–307. https://doi.org/10.1590/S0102-30982014000200004.

Siegel, J. S., and D. A. Swanson. 2004. The Methods and Materials of Demography. Bingley, England, UK: Emerald Group Publishing Limited. https://books.google.com.ar/books/about/The_Methods_and_Materials_of_Demography.html?id=eJDG1OeVLxoC&redir_esc=y.

United Nations, Department of Economic, and Population Division Social Affairs. 2022. “World Population Prospects 2022: Methodology of the United Nations population estimates and projections” DESA/POP/2022/TR/NO. 4. https://population.un.org/wpp/Publications/Files/WPP2022_Methodology.pdf.

Wachter, Kenneth W. 2006. “Essential Demographic Methods.”

Wilmoth, John, Sarah Zureick, Vladimir Canudas-Romo, Mie Inoue, and Cheryl Sawyer. 2012. “A Flexible Two-Dimensional Mortality Model for Use in Indirect Estimation.” Population Studies 66 (1): 1. https://doi.org/10.1080/00324728.2011.611411.

UNDESA-División de Población, ivan.williams@un.org ↩︎

Análisis demográfico en R: técnicas básicas y DemoTools

Estimaciones de población y ronda censal 2020: Desafíos y lecciones aprendidas en tiempos de pandemia

Iván Williams¹

Actualizado: 2023-11-23

Introducción

Objetivo y setting inicial

Data

Población

Evolución censal

Edad y sexo

Edad ignorada

Atracción de edad y suavizamiento

Edades iniciales

Mortalidad

Construcción y graduación

Uso de tablas modelo

Fecundidad

Graduación

Otros

Procedimiento de ajuste proporcional iterativo

Patrones de migración neta por edad

Actividad

Referencias

Análisis demográfico en R: técnicas básicas y DemoTools

Estimaciones de población y ronda censal 2020: Desafíos y lecciones aprendidas en tiempos de pandemia

Iván Williams1

Actualizado: 2023-11-23

Introducción

Objetivo y setting inicial

Data

Población

Evolución censal

Edad y sexo

Edad ignorada

Atracción de edad y suavizamiento

Edades iniciales

Mortalidad

Construcción y graduación

Uso de tablas modelo

Fecundidad

Graduación

Otros

Procedimiento de ajuste proporcional iterativo

Patrones de migración neta por edad

Actividad

Referencias

Iván Williams¹