21  Стандартизированные коэффициенты

На этой странице мы покажем два способа стандартизировать исход, например, госпитализации или смертность, по таким характеристикам как возраст и пол.

• Используя пакет dsr
• Используя пакет PHEindicatormethods

Мы начнем с подробной демонстрации процесса подготовки/вычистки/соединения данных, так как это часто требуется при объединении популяционных данных из разных стран, данных о стандартном населении, смертях и т.п.

21.1 Обзор

Существует два основных способа стандартиации: прямая и косвенная стандартизация Представим, что вам нужно стандартизировать коэффициент смертности по возрасту и полу для страны A и страны B, а также сравнить стандартизированные коэффициенты между этими странами.

• Для прямой стандартизации необходимо знать численность населения, подверженного риску, и количество смертей для каждой половозрастной страты в стране А и стране Б. Одной из страт в нашем примере могут быть женщины в возрасте 15-44 лет.
  
• Для косвенной стандартизации достаточно знать общее число умерших и половозрастную структуру населения каждой страны. Поэтому этот вариант возможен, если отсутствуют возрастные и половые коэффициенты смертности или численность населения. Кроме того, косвенная стандартизация предпочтительна в случае небольшого числа страт, так как на оценки при прямой стандартизации будет влиять значительный разброс выборки.

21.2 Подготовка

Чтобы показать, как проводится стандартизация, мы будем использовать фиктивные данные о численности населения и смертности в странах А и Б с разбивкой по возрасту (в категориях 5 лет) и полу (женщины, мужчины). Чтобы подготовить наборы данных к использованию, выполним следующие подготовительные действия:

1. Загрузим пакеты
   
2. Загрузим наборы данных
   
3. Соединим данные по населению и смертям по двум странам
4. Повернем вертикально, чтобы была одна строка на половозрастную страту
5. Вычистим референтную популяцию (мировое стандартное население) и присоединим ее к страновым данным

В вашем сценарии данные могут быть представлены в другом формате. Возможно, данные представлены по провинциям, городам или другим районам. У вас может быть одна строка для каждого случая смерти и информация о возрасте и поле для каждого (или значительной части) из этих случаев смерти. В таком случае, см. страницы Группирование данных, Поворот данных и Описательные таблицы, чтобы создать набор данных с подсчетом событий и населения по половозрастным стратам.

Нам также нужна референтная популяция, стандартное население. В этом упражнении мы будет использовать world_standard_population_by_sex. Мировое стандартное население основано на населении 46 стран и было разработано в 1960. Существует много “стандартных” популяций - в качестве примера см. веб-сайт Медицинской службы Шотландии, который является очень информативным в вопросах Европейского стандартного населения, Мирового стандартного населения и Стандартного населения Шотландии.

Загрузка пакетов

Данный фрагмент кода показывает загрузку пакетов, необходимых для анализа. В данном руководстве мы фокусируемся на использовании p_load() из пакета pacman, которая устанавливает пакет, если необходимо, и загружает его для использования. Вы можете также загрузить установленные пакеты с помощью library() из базового R. См. страницу Основы R для получения дополнительной информации о пакетах R.

pacman::p_load(
     rio,                 # импорт/экспорт данных
     here,                # путь к файлам
     stringr,             # вычистка текста и последовательностей
     frailtypack,         # необходим для dsr, для моделей, учитывающих индивидуальную уязвимость
     dsr,                 # стандартизированные коэффициенты
     PHEindicatormethods, # альтернатива для стандартизации коэффициентов
     tidyverse)           # управление данными и визуализация

ВНИМАНИЕ: Если у вас более новая версия R, пакет dsr нельзя напрямую скачать с CRAN. Однако он все еще доступен в архиве CRAN. Вы можете установить и использовать этот.

Дляф всех пользователей, кроме Mac:

packageurl <- "https://cran.r-project.org/src/contrib/Archive/dsr/dsr_0.2.2.tar.gz"
install.packages(packageurl, repos=NULL, type="source")

# другое решение, которое может сработать
require(devtools)
devtools::install_version("dsr", version="0.2.2", repos="http:/cran.us.r.project.org")

Для пользователей Mac:

require(devtools)
devtools::install_version("dsr", version="0.2.2", repos="https://mac.R-project.org")

Загружаем популяционные данные

Для получения инструкций о том, как скачать все данные для примеров из руководства см. страницу [Скачивание руководства и данных]. Вы можете импортировать данные для страницы Стандартизация напрямую в R из нашего репозитория Github, выполнив следующие команды import():

# импорт демографических данных по стране A напрямую из Github
A_demo <- import("https://github.com/appliedepi/epirhandbook_eng/raw/master/data/standardization/country_demographics.csv")

# импорт данных о смертях для страны A напрямую из Github
A_deaths <- import("https://github.com/appliedepi/epirhandbook_eng/raw/master/data/standardization/deaths_countryA.csv")

# импорт демографических данных по стране B напрямую из Github
B_demo <- import("https://github.com/appliedepi/epirhandbook_eng/raw/master/data/standardization/country_demographics_2.csv")

# импорт данных о смертях для страны B напрямую из Github
B_deaths <- import("https://github.com/appliedepi/epirhandbook_eng/raw/master/data/standardization/deaths_countryB.csv")

# импорт демографических данных по стране B напрямую из Github
standard_pop_data <- import("https://github.com/appliedepi/epirhandbook_eng/raw/master/data/standardization/world_standard_population_by_sex.csv")

Сначала мы загружаем демографические данные (количество мужчин и женщин по 5-летним возрастным категориям) для двух стран, которые мы будем сравнивать, “Страны A” и “Страны B”.

# Страна A
A_demo <- import("country_demographics.csv")

# Страна B
B_demo <- import("country_demographics_2.csv")

Загружаем данные о количестве смертей

Что удобно, у нас также естьт количество смертей за интересующий временной период по полу и возрасту. Количество по каждой стране представлено в отдельном файле, показанном ниже.

Смерти в Стране A

Смерти в Стране B

Вычищаем данные по населению и смертям

Нам нужно соединить и трансформировать эти данные следующим образом:

• Объединить населения стран в один набор данных и повернуть вертикально, чтобы каждая половозрастная страта занимала одну строку
  
• Объединить количество смертей стран и повернуть вертикально, чтобы каждая половозрастная страта занимала одну строку
  
• Присоединить смерти к популяциям

Сначала объединим наборы данных с населением стран, повернем вертикально и проведем небольшую вычистку. См. страницу [Поворот данных] для получения более подробной информации.

pop_countries <- A_demo %>%  # начнем с набора данных по Стране A
     bind_rows(B_demo) %>%        # связываем строки, поскольку столбцы имеют одинаковые названия
     pivot_longer(                       # поворачиваем вертикально
          cols = c(m, f),                   # столбцы объединяются в один
          names_to = "Sex",                 # имя для нового столбца, содержащего категорию ("m" или "f") 
          values_to = "Population") %>%     # имя для нового столбца, содержащего повернутые числовые значения
     mutate(Sex = recode(Sex,            # перекодируем значения для ясности
          "m" = "Male",
          "f" = "Female"))

Объединенные данные по населению теперь выглядят вот так (пролистайте, чтобы увидеть страны A и B):

Теперь проведем похожие операции с двумя наборами данных о смертях.

deaths_countries <- A_deaths %>%    # начинаем с набора данных о смертях в стране A
     bind_rows(B_deaths) %>%        # связываем строки с набором данных B, поскольку столбцы имеют одинаковые названия
     pivot_longer(                  # поворачиваем вертикально
          cols = c(Male, Female),        # столбцы трансформируются в один
          names_to = "Sex",              # имя для нового столбца, содержащего категорию ("m" или "f") 
          values_to = "Deaths") %>%      # имя для нового столбца, содержащего повернутые числовые значения
     rename(age_cat5 = AgeCat)      # переименуем для ясности

Данные о смертях теперь выглядят вот так и содержат данные из обеих стран:

Теперь мы соединим данные о смертях и о населении на основе общих столбцов Country (страна), age_cat5 (возрастная категория) и Sex (пол). Это добавит столбец Deaths (смерти).

country_data <- pop_countries %>% 
     left_join(deaths_countries, by = c("Country", "age_cat5", "Sex"))

Теперь мы классифицируем Sex, age_cat5 и Country как факторы и установим порядок уровней с помощью функции fct_relevel() из пакета forcats, как описано на странице [Факторы]. Обратите внимание, что классификация уровней фактора не меняет видимые данные, а команда arrange() сортирует по Стране, возрастной категории и полу.

country_data <- country_data %>% 
  mutate(
    Country = fct_relevel(Country, "A", "B"),
      
    Sex = fct_relevel(Sex, "Male", "Female"),
        
    age_cat5 = fct_relevel(
      age_cat5,
      "0-4", "5-9", "10-14", "15-19",
      "20-24", "25-29",  "30-34", "35-39",
      "40-44", "45-49", "50-54", "55-59",
      "60-64", "65-69", "70-74",
      "75-79", "80-84", "85")) %>% 
          
  arrange(Country, age_cat5, Sex)

ВНИМАНИЕ: Если у вас мало смертей на страту, рассмотрите возможность использования 10- или 15-летних категорий вместо 5-летних категорий возраста.

Загрузка референтной популяции

И последнее, для прямой стандартизации мы импортируем референтную популяцию (мировое “стандартное население” по полу)

# референтная популяция
standard_pop_data <- import("world_standard_population_by_sex.csv")

Вычистка референтной популяции

Значения возрастных категорий в country_data и standard_pop_data необходимо унифицировать.

В настоящее время значения столбца age_cat5 в датафрейме standard_pop_data содержат слова “years” и “plus”, а в датафрейме country_data их нет. Нам нужно, чтобы значения возрастных категорий совпадали. Мы используем str_replace_all() из пакета stringr, как описано на странице [Текст и последовательности], чтобы заменить эти комбинации символов на отсутствие пробела "".

Более того, пакет dsr ожидает, что в стандартном населении столбец, где содержится количество, должен называться "pop". Поэтому нам нужно переименовать столбец соответствующим образом.

# удаляем конкретную последовательность из значений столбца
standard_pop_clean <- standard_pop_data %>%
     mutate(
          age_cat5 = str_replace_all(age_cat5, "years", ""),   # удаляем "year"
          age_cat5 = str_replace_all(age_cat5, "plus", ""),    # удаляем "plus"
          age_cat5 = str_replace_all(age_cat5, " ", "")) %>%   # удаляем пробел " "
     
     rename(pop = WorldStandardPopulation)   # меняем имя столбца на "pop", так как этого ожидает пакет dsr

ВНИМАНИЕ: Если вы будете пытаться использовать str_replace_all(), чтобы удалить символ плюс, она не сработает, поскольку это специальный символ. “Изолируйте” этот особый характер символа, поставив перед ним два обратных слэша, вот так str_replace_call(column, "\\+", "").

Создание набора данных со стандартным населением

Наконец пакет PHEindicatormethods, который де6тально описан ниже, ожидает присоединения стандартных популяций к данным о количестве случаев события и населении. Поэтому в этих целях мы создадим набор данных all_data.

all_data <- left_join(country_data, standard_pop_clean, by=c("age_cat5", "Sex"))

Этот полный набор данных выглядит следующим образом:

21.3 Пакет dsr

Ниже мы демонстрируем расчет и сравнение прямых стандартизированных коэффициентов, используя пакет dsr. Пакет dsr позволяет вам рассчитывать и сравнивать прямые стандартизированные коэффициенты (не косвенные стандартизированные коэффициенты!).

В разделе Подготовка данных мы создали отдельные наборы данных для подсчета по странам и стандартного населения:

1. объект country_data, который является популяционной таблицей с количеством населения и количеством смертей по стратам по стране
   
2. объект standard_pop_clean, содержащий количество населения по стратам для референтной популяции, Мирового стандартного населения

Мы будем использовать эти отдельные наборы данных для подхода с использованием dsr.

Стандартизированные коэффициенты

Ниже мы рассчитываем коэффициенты по стране, стандартизированные прямым образом по полу и возрасту. Мы используем функцию dsr().

Необходимо отметить, что dsr() ожидает одного датафреймя для населения страны и количество случаев события (смертей) и отдельный датафрейм с референтной популяцией. Она также ожидает, что в этом наборе референтной популяции имя столбца за единицу времени будет “pop” (мы это сделали на этапе Подготовки данных).

Существует много аргументов, как показано в коде ниже. Необходимо отметить, что event = задается как столбец Deaths, а fu = (“отслеживание”) задается как столбец Population. Мы устанавливаем подгруппы сравнения как столбец Country и стандартизируем на основе age_cat5 и Sex. Этим последним двум столбцам не присваивается конкретный именованный аргумент. См. ?dsr для получения дополнительной информации.

# Рассчитываем коэффициенты по стране, стандартизированные прямым образом по возрасту и полу
mortality_rate <- dsr::dsr(
     data = country_data,  # уточняем объект, содержащий количество смертей по страте
     event = Deaths,       # столбец, содержащий количество смертей по страте 
     fu = Population,      # столбец, содержащий количество населения по страте
     subgroup = Country,   # единицы, которые мы бы хотели сравнить
     age_cat5,             # другие столбцы - коэффициенты будут стандартизированы по этим
     Sex,
     refdata = standard_pop_clean, # датафрейм референтной популяции со столбцом под названием pop
     method = "gamma",      # метод расчета 95% ДИ
     sig = 0.95,            # уровень значимости
     mp = 100000,           # нам нужны коэффициенты на 100 000 населения
     decimals = 2)          # количество десятичных знаков)


# Печать выходнрых данных в виде красивой таблицы HTML
knitr::kable(mortality_rate) # показываем коэффициент смертности до и после прямой стандартизации

Subgroup	Numerator	Denominator	Crude Rate (per 1e+05)	95% LCL (Crude)	95% UCL (Crude)	Std Rate (per 1e+05)	95% LCL (Std)	95% UCL (Std)
A	11344	86790567	13.07	12.83	13.31	23.57	23.08	24.06
B	9955	52898281	18.82	18.45	19.19	19.33	18.46	20.22

Выше мы видим, что хотя у страны A был ниже общий коэффициент смертности, чем у страны B, у нее выше стандартизированный коэффициент после прямой стандартизации по полу и возрасту.

Standardized rate ratios

# Рассчитываем ОР
mortality_rr <- dsr::dsrr(
     data = country_data, # уточняем объект, содержащий количество смертей по страте
     event = Deaths,      # столбец, содержащий количество смертей по страте 
     fu = Population,     # столбец, содержащий количество населения по страте
     subgroup = Country,  # единицы, которые мы бы хотели сравнить
     age_cat5,
     Sex,                 # характеристики, по которым мы хотим стандартизировать 
     refdata = standard_pop_clean, # референтная популяция, с числами в столбце под названием pop
     refgroup = "B",      # референс для сравнения
     estimate = "ratio",  # тип оценки
     sig = 0.95,          # уровень значимости
     mp = 100000,         # нам нужны коэффициенты на 100 000 населения
     decimals = 2)        # количество десятичных знаков

# Печать таблицы
knitr::kable(mortality_rr) 

Comparator	Reference	Std Rate (per 1e+05)	Rate Ratio (RR)	95% LCL (RR)	95% UCL (RR)
A	B	23.57	1.22	1.17	1.27
B	B	19.33	1.00	0.94	1.06

Стандартизированный коэффициент смертности в 1.22 раза выше в стране A по сравнению со страной B (95% ДИ 1.17-1.27).

Разница стандартизированных коэффициентов

# Calculate RD
mortality_rd <- dsr::dsrr(
     data = country_data,       # уточняем объект, содержащий количество смертей по страте
     event = Deaths,            # столбец, содержащий количество смертей по страте 
     fu = Population,           # столбец, содержащий количество населения по страте
     subgroup = Country,        # единицы, которые мы бы хотели сравнить
     age_cat5,                  # характеристики, по которым мы хотим стандартизировать
     Sex,                        
     refdata = standard_pop_clean, # референтная популяция, с числами в столбце под названием pop
     refgroup = "B",            # референс для сравнения
     estimate = "difference",   # тип оценки
     sig = 0.95,                # уровень значимости
     mp = 100000,               # нам нужны коэффициенты на 100 000 населения
     decimals = 2)              # количество десятичных знаков

# Печать таблицы
knitr::kable(mortality_rd) 

Comparator	Reference	Std Rate (per 1e+05)	Rate Difference (RD)	95% LCL (RD)	95% UCL (RD)
A	B	23.57	4.24	3.24	5.24
B	B	19.33	0.00	-1.24	1.24

В стране A имеется 4.24 дополнительных смертей на 100.000 населения (95% ДИ 3.24-5.24) по сравнению со страной A.

21.4 Пакет PHEindicatormethods

Еще один способ расчета стандартизированных коэффициентов - использовать пакет PHEindicatormethods. Этот пакет позволяет вам рассчитывать как прямые, так и косвенные стандартизированные коэффициенты. Мы покажем оба варианта.

В этом разделе используется датафрейм all_data, созданный в конце раздела Подготовка. Этот датафрейм включает населения стран, смерти и мировую стандартную референтную популяцию. Вы можете посмотреть это тут.

Прямые стандартизированные коэффициенты

Ниже мы сначала группируем данные по Стране, затем передаем их в функцию phe_dsr(), чтобы получить прямые стандартизированные коэффициенты по странам.

Следует отметить, что референтное (стандартное) население может быть указано как столбец внутри датафрейма для конкретной страны, либо как отдельный вектор. Если вы его указываете внутри странового датафрейма, вы должны установить аргумент stdpoptype = "field". Если вы его указываете как вектор, установите stdpoptype = "vector". В последнем случае вам необходимо убедиться, что упорядочивание строк по стратам аналогично и в страновом датафрейме, и в референтной популяции, так как записи будут сопоставляться по позиции. В примере ниже мы задали референтную популяцию как столбец внутри странового датафрейма.

См. справку с помощью ?phr_dsr или ссылки в разделе Ресурсы для получения более подробной информации.

# Рассчитываем коэффициенты по стране, стандартизированные прямым образом по полу и возрасту
mortality_ds_rate_phe <- all_data %>%
     group_by(Country) %>%
     PHEindicatormethods::phe_dsr(
          x = Deaths,                 # столбец с наблюдаемым количеством событий
          n = Population,             # столбец с нестандартным населением для каждой страты
          stdpop = pop,               # стандартное население для каждой страты
          stdpoptype = "field")       # либо "vector" (вектор) для отдельного вектора, либо "field" (поле), что означает, что стандартное население указано в данных  

# Печать таблицы
knitr::kable(mortality_ds_rate_phe)

Country	total_count	total_pop	value	lowercl	uppercl	confidence	statistic	method
A	11344	86790567	23.56686	23.08107	24.05944	95%	dsr per 100000	Dobson
B	9955	52898281	19.32549	18.45516	20.20882	95%	dsr per 100000	Dobson

Косвенные стандартизированные коэффиенты

Для косвенной стандартизации вам нужна референтная популяция с количеством смертей и количеством населения на страту. В данном примере мы будем рассчитывать коэффициенты для Страны A, используя страну B в качестве референтной популяции, так как референтная популяция standard_pop_clean не включает количество смертей на страту.

Ниже мы сначала создаем референтную популяцию из страны B. Затем мы передаем данные о смертности и населении страны A, объединяем их с референтной популяцией и передаем в функцию calculate_ISRate(), чтобы получить косвенные стандартизированные коэффициенты. Конечно, можно сделать и наоборот.

Следует отметить - в примере ниже референтная популяция задается как отдельный датафрейм. В этом случае нужно убедиться, что векторы x =, n =, x_ref = и n_ref = все упорядочены по тем же значениям категории стандартизации (страты), что и в страновом датафрейме, так как сопоставление записей делается по позиции.

См. справку ?phr_isr или ссылки в разделе Ресурсы для получения дополнительной информации.

# Создаем референтную популяцию
refpopCountryB <- country_data %>% 
  filter(Country == "B") 

# Рассчитываем коэффициенты для страны А, косвенным образом стандартизированные по возрасту и полу
mortality_is_rate_phe_A <- country_data %>%
     filter(Country == "A") %>%
     PHEindicatormethods::calculate_ISRate(
          x = Deaths,                 # столбец с наблюдаемым количеством событий
          n = Population,             # столбец с нестандартным населением для каждой страты
          x_ref = refpopCountryB$Deaths,  # референтное количество смертей для каждой страты
          n_ref = refpopCountryB$Population)  # референтная популяция по каждой страте

# Печать таблицы
knitr::kable(mortality_is_rate_phe_A)

observed	expected	ref_rate	value	lowercl	uppercl	confidence	statistic	method
11344	15847.42	18.81914	13.47123	13.22446	13.72145	95%	indirectly standardised rate per 100000	Byars

21.5 Ресурсы

Если вам нужен еще один воспроизводимый пример использования dsr, см. эту виньетку

Еще один пример использования PHEindicatormethods доступен на этом сайте

См. PHEindicatormethods справочный pdf файл