Демонстрация навыков

Author

Погожий Ю.С.

О ЧЕМ РЕЧЬ?

Продемонстрированы этапы статистической обработки и анализа информации историй болезни пациентов, госпитализированных в инфекционный центр в период с апреля 2020 г. по февраль 2021 г. с помощью программного обеспечения R-studio (язык программирования R).

ПЕРВЫЙ ЭТАП: ИМПОРТИРОВАНИЕ БАЗЫ ДАННЫХ ЭЛЕКТРОННЫХ ЗАПИСЕЙ ИСТОРИЙ БОЛЕЗНЕЙ.

### Загрузка необходимых расширений
options(scipen = 999) 
x<-c("tidyverse", "lubridate", "readxl", "openxlsx", "finalfit",
     "gapminder", "broom")
sapply(x, require, character.only=TRUE)

tidyverse lubridate    readxl  openxlsx  finalfit gapminder     broom 
     TRUE      TRUE      TRUE      TRUE      TRUE      TRUE      TRUE

### Импортирование таблицы Exel, первичная визуализация

data <- read_xlsx("/Users/insomnia/R/статья/data.xlsx")

### Очистка и форматирование табличных данных

time <- 2020
dataset <- data %>% mutate(Birthdate = year(Birthdate)) %>% mutate(age = time - Birthdate)
dataset <- dataset %>% mutate(discharge = ymd(discharge)) 
dataset <- dataset %>% mutate(contingent = as.factor(contingent), outcome = as.factor(outcome), `Disease severity` = as.factor(`Disease severity`), 
                              `respiratory support` = as.factor(`respiratory support`), `ICU admission` = as.factor(`ICU admission`), Diagnosis = as.factor(Diagnosis))
dataset <- dataset %>% select(-discharge)
dataset <- dataset %>% select(Name, historyN, Birthdate, age, contingent, Diagnosis, everything()) %>% na.exclude()
dataset <- dataset %>% mutate(`Disease severity` = fct_relevel(.f = dataset$`Disease severity`, "mild", "moderate", "severe", "critical"))

### Исключение неинформативных столбцов

dataset_clean <- dataset %>% select(-Name, -historyN, -Birthdate, -contingent, -Diagnosis)
colnames <- c("age", "hosp_days", "outcome", "Disease_severity", "respiratory_support", "ICU_admission")
colnames(dataset_clean) <- colnames
dataset_clean <- dataset_clean %>% filter(outcome %in% c("умер", "выбыл с улучшением"))
dataset_clean$outcome <- factor(dataset_clean$outcome, levels = c("умер", "выбыл с улучшением"))
### Результат. Рабочая таблица

print(dataset_clean)

# A tibble: 257 × 6
     age hosp_days outcome    Disease_severity respiratory_support ICU_admission
   <dbl>     <dbl> <fct>      <fct>            <fct>               <fct>        
 1    60        12 выбыл с у… mild             none                no           
 2    56        16 выбыл с у… severe           non-invasive        yes          
 3    51        16 выбыл с у… moderate         none                no           
 4    72        10 выбыл с у… moderate         non-invasive        no           
 5    47        22 выбыл с у… moderate         none                no           
 6    59        15 выбыл с у… severe           non-invasive        yes          
 7    48        13 выбыл с у… severe           non-invasive        yes          
 8    48        16 выбыл с у… severe           non-invasive        yes          
 9    40        12 выбыл с у… moderate         non-invasive        no           
10    78        17 выбыл с у… moderate         non-invasive        yes          
# ℹ 247 more rows

ВЛИЯЕТ ЛИ ВОЗРАСТ НА СМЕРТНОСТЬ?

### Общая информация по группам

dataset_clean %>% group_by(outcome) %>% summarise("возраст" = mean(age), n = n()) %>% 
  mutate("смертность" = round(n / sum(n), 3))

# A tibble: 2 × 4
  outcome            возраст     n смертность
  <fct>                <dbl> <int>      <dbl>
1 умер                  80.1    17      0.066
2 выбыл с улучшением    60.9   240      0.934

### тест Шапиро-Уилка. Используется для определения того, соответствует ли выборка нормальному распределению.

shapiro.test(dataset_clean$age)


    Shapiro-Wilk normality test

data:  dataset_clean$age
W = 0.98151, p-value = 0.002013

### Значение p теста оказывается равным 0,002013. Поскольку это значение меньше 0,05, у нас есть достаточно доказательств, чтобы сказать, что данные выборки не получены из населения с нормальным распределением.

### Визуализируем данные

dataset_clean %>% ggplot(aes(x = outcome, y = age, fill = outcome)) +
  geom_boxplot()

dataset_clean %>%
  ggplot(aes(x = age)) +
  geom_density() +
  facet_wrap(~ outcome)

dataset_clean %>%
  ggplot(aes(x = age,
             y = outcome,
             fill = outcome)) +
  ggridges::geom_density_ridges(bandwidth = 6)

### U-критерий Манна — Уитни — Уилкоксона (англ. Mann–Whitney–Wilcoxon), — статистический критерий, используемый для оценки различий между двумя независимыми выборками по уровню какого-либо признака, измеренного количественно. Позволяет выявлять различия в значении параметра между малыми выборками.

dataset_clean %>% rstatix::wilcox_test(age ~ outcome, paired = FALSE)

# A tibble: 1 × 7
  .y.   group1 group2                n1    n2 statistic          p
* <chr> <chr>  <chr>              <int> <int>     <dbl>      <dbl>
1 age   умер   выбыл с улучшением    17   240     3472. 0.00000132

### Различия есть, p = 0.00000132

ЗАВИСИТ ЛИ СТЕПЕНЬ ТЯЖЕСТИ ЗАБОЛЕВАНИЯ ОТ ВОЗРАСТА?

### 4 группы:

dataset_clean %>% group_by(Disease_severity) %>% summarise(mean = mean(age), n = n())

# A tibble: 4 × 3
  Disease_severity  mean     n
  <fct>            <dbl> <int>
1 mild              56.6    78
2 moderate          60.9    83
3 severe            65.8    60
4 critical          71.0    36

lab_names <- list(
  'mild'="легкая",
  'moderate'="средняя",
  'severe'="тяжелая",
  'critical'="кр. тяжелая"
)
hospital_labeller <- function(variable,value){
  return(lab_names[value])
}

### Гистограмма: возраст у 4 групп пациентов в зависимости от тяжести заболевания

dataset_clean %>% ggplot(aes(x = age)) +
  geom_histogram(bins = 20) +
  facet_grid(~Disease_severity, labeller = hospital_labeller)+
  labs(x = "возраст",
       y = "количество",
       title = "Гистограмма")+
  theme(plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5))

### График квантиль-квантиль для визуального выявления отклонений от нормальности

dataset_clean %>% ggplot(aes(sample = age)) +
  geom_qq() +
  geom_qq_line(colour = "blue") +
  facet_grid(~Disease_severity, labeller = hospital_labeller) +
labs(title = "График квантиль-квантиль")+
  theme(plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5))

СРАВНЕНИЕ КОЛИЧЕСТВЕННОГО ПРИЗНАКА (ВОЗРАСТ) ВНУТРИ ГРУПП

т.к. групп 4 (более 2) и распределение отличается от нормального, используется непараметрический критерий Краскела-Уоллиса.

kruskal.test(age ~ Disease_severity, data = dataset_clean)


    Kruskal-Wallis rank sum test

data:  age by Disease_severity
Kruskal-Wallis chi-squared = 22.849, df = 3, p-value = 0.00004343

### p-value = 0,00004343, различия есть, но между какими группами?

### В случае выявления различий при сравнении всех групп с помощью критерия Краскела-Уоллиса для более точного описания наблюдаемых тенденций используется критерии, позволяющие оценить различия показателей при сравнении групп попарно. Для этого подойдет парный критерий Манна-Уитни-Уилкоксона.

pairwise.wilcox.test(dataset_clean$age, dataset_clean$Disease_severity,
                 p.adjust.method = "BH")


    Pairwise comparisons using Wilcoxon rank sum test with continuity correction 

data:  dataset_clean$age and dataset_clean$Disease_severity 

         mild     moderate severe
moderate 0.1355   -        -     
severe   0.0085   0.1355   -     
critical 0.000078 0.0015   0.0648

P value adjustment method: BH

### Диаграмма размаха («ящик с усами») (Box and Whisker Plot или Box Plot)

dataset_clean %>% ggplot(aes(y = age, x = `Disease_severity`)) +
  geom_boxplot()+
  scale_x_discrete(labels = c("легкая", "средняя", "тяжелая", "крайне тяжелая"))+
  labs(y = "возраст",
       x = "тяжесть заболевания", 
       title = "Box Plot") +
  theme(plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5))

ВЫВОДЫ:

1. С увеличением возраста возрастает смертность

2. При парном сравнении групп (‘mild’=“легкая”, ‘moderate’=“средняя”, ‘severe’=“тяжелая”, ‘critical’=“кр. тяжелая”), статистически значимые различия (p<0.05) обнаружены между группами:

легкая - тяжелая, легкая - кр. тяжелая, средняя - крайне тяжелая.

Статистически незначимые различия между группами:

легкая - средняя, средняя - тяжелая, тяжелая - кр. тяжелая.

broom::tidy(pairwise.wilcox.test(dataset_clean$age, dataset_clean$Disease_severity,
                 p.adjust.method = "BH"))

# A tibble: 6 × 3
  group1   group2     p.value
  <chr>    <chr>        <dbl>
1 moderate mild     0.136    
2 severe   mild     0.00846  
3 severe   moderate 0.136    
4 critical mild     0.0000785
5 critical moderate 0.00153  
6 critical severe   0.0648