Gọi package

Đoạn code này hiển thị những package cần tải cho các phân tích. Trong sổ tay này, chúng tôi nhấn mạnh đến hàm p_load() từ pacman, hàm sẽ cài đặt package nếu cần và gọi nó ra để sử dụng. Bạn cũng có thể gọi các package đã cài đặt với library() từ base R. Xem chương R cơ bản để có thêm thông tin về các R package.

pacman::p_load(
  tidyverse,      # for data management and viz
  slider,         # for calculating moving averages
  tidyquant       # for calculating moving averages within ggplot
)

Nhập dữ liệu

Chúng ta nhập các trường hợp trong linelist đã được làm sạch từ một vụ dịch Ebola mô phỏng. Nếu bạn muốn theo dõi, bấm để tải xuống linelist “đã được làm sạch” (tệp .rds). Nhập dữ liệu với hàm import() từ package rio (hàm này xử lý nhiều loại tệp như .xlsx, .csv, .rds - xem chương Nhập xuất dữ liệu để biết thêm chi tiết).

# import the linelist
linelist <- import("linelist_cleaned.xlsx")

50 hàng đầu tiên của linelist được hiển thị dưới đây.

Cuộn theo ngày

Trường hợp sử dụng có khả năng xảy ra nhất của tính toán biến động trong dịch tễ học ứng dụng là kiểm tra một số liệu theo thời gian. Ví dụ: đo lường động các ca mới mắc, dựa trên số lượng trường hợp hàng ngày.

Nếu bạn có dữ liệu chuỗi thời gian đã được làm sạch với đủ giá trị cho tất cả các ngày, bạn có thể sử dụng hàm slide_dbl(), như đã được trình bày trong chương Chuỗi thời gian và phát hiện ổ dịch.

Tuy nhiên, trong nhiều trường hợp dịch tễ học ứng dụng, bạn có thể gặp những ngày trống trong dữ liệu của mình, những ngày mà không có sự kiện nào được ghi lại. Trong những trường hợp này, tốt nhất là sử dụng các phiên bản “chỉ mục” của các hàm slider.

Dữ liệu được lập chỉ mục

Dưới đây, chúng tôi trình bày một ví dụ sử dụng slide_index_dbl() đối với các trường hợp của bộ dữ liệu linelist. Giả sử rằng mục tiêu của chúng ta là tính toán tỷ lệ mới mắc liên tục trong 7 ngày - tính tổng các trường hợp bằng cách sử dụng cửa sổ 7 ngày luân phiên. Nếu bạn đang tìm kiếm ví dụ về trung bình động, hãy xem phần bên dưới về cuộn theo nhóm.

Để bắt đầu, bộ dữ liệu daily_counts được tạo ra để phản ánh số lượng ca mắc hàng ngày từ linelist, như đã được tính toán với hàm count() trong dplyr.

# make dataset of daily counts
daily_counts <- linelist %>% 
  count(date_hospitalisation, name = "new_cases")

Dưới đây là data frame daily_counts - bao gồm nrow(daily_counts) hàng, mỗi hàng đại diện cho một ngày, nhưng đặc biệt trong giai đoạn đầu của dịch, có một số ngày không xuất hiện (không có ca mắc nào được tiếp nhận vào những ngày đó).

Điều quan trọng là phải nhận ra rằng một hàm cuộn tiêu chuẩn (như slide_dbl() sẽ sử dụng cửa sổ của 7 hàng, không phải 7 ngày. Vì vậy, nếu có bất kỳ ngày nào trống, một số cửa sổ sẽ thực sự kéo dài hơn 7 ngày theo lịch!

Một cửa sổ động “thông minh” có thể được tạo với hàm slide_index_dbl(). “Chỉ mục” có nghĩa là hàm sử dụng một cột riêng biệt làm “chỉ mục” cho cửa sổ động. Cửa sổ đó không chỉ đơn giản dựa trên các hàng của data frame.

Nếu cột chỉ mục là ngày, bạn có thêm khả năng xác định phạm vi cửa sổ cho .before = và/hoặc .after = theo đơn vị days() và months() của lubridate. Nếu bạn thực hiện những điều này, hàm sẽ bao gồm những ngày trống trong cửa sổ như thể chúng ở đó (dưới dạng giá trị NA).

Hãy đưa ra một so sánh. Dưới đây, chúng tôi tính toán số trường hợp mới mắc biến động trong 7 ngày với các cửa sổ thông thường và được lập chỉ mục.

rolling <- daily_counts %>% 
  mutate(                                # create new columns
    # Using slide_dbl()
    ###################
    reg_7day = slide_dbl(
      new_cases,                         # calculate on new_cases
      .f = ~sum(.x, na.rm = T),          # function is sum() with missing values removed
      .before = 6),                      # window is the ROW and 6 prior ROWS
    
    # Using slide_index_dbl()
    #########################
    indexed_7day = slide_index_dbl(
        new_cases,                       # calculate on new_cases
        .i = date_hospitalisation,       # indexed with date_onset 
        .f = ~sum(.x, na.rm = TRUE),     # function is sum() with missing values removed
        .before = days(6))               # window is the DAY and 6 prior DAYS
    )

Quan sát trong cột thông thường cho 7 hàng đầu tiên, cách số lượng ca mắc tăng đều đặn mặc dù các hàng không cách nhau 7 ngày! Cột “được lập chỉ mục” liền kề tính cho những ngày lịch trống này, vì vậy tổng 7 ngày của nó thấp hơn nhiều, ít nhất là khi khoảng thời gian các ca mắc cách nhau xa hơn trong thời kỳ dịch bệnh này.

Bây giờ bạn có thể vẽ biểu đồ những dữ liệu này bằng cách sử dụng ggplot():

ggplot(data = rolling)+
  geom_line(mapping = aes(x = date_hospitalisation, y = indexed_7day), size = 1)

Biến động theo nhóm

Nếu bạn nhóm dữ liệu của mình trước khi sử dụng hàm slider, các cửa sổ trượt sẽ được áp dụng theo nhóm. Hãy cẩn thận để sắp xếp các hàng của bạn với thứ tự mong muốn theo nhóm.

Mỗi khi một nhóm mới được tạo, cửa sổ trượt sẽ bắt đầu lại. Do đó, một điều cần lưu ý là nếu dữ liệu của bạn được nhóm lại và bạn đã thiết lập .complete = TRUE, bạn sẽ có các giá trị trống ở mỗi lần dịch chuyển giữa các nhóm. Khi hàm di chuyển xuống dưới qua các hàng, mọi dịch chuyển trong cột được nhóm sẽ bắt đầu lại việc cộng dồn kích thước cửa sổ tối thiểu để cho phép tính toán.

Xem chương Nhóm dữ liệu trong sổ tay này để biết thêm chi tiết về nhóm dữ liệu.

Dưới đây, chúng tôi đếm các ca mắc trong linelist theo ngày và theo bệnh viện. Sau đó, chúng tôi sắp xếp các hàng theo thứ tự tăng dần, thứ tự đầu tiên theo bệnh viện và sau đó là theo ngày. Tiếp theo, chúng tôi đặt group_by(). Cuối cùng, chúng tôi có thể tạo trung bình động mới của mình.

grouped_roll <- linelist %>%
     
  count(hospital, date_hospitalisation, name = "new_cases") %>% 
     
  arrange(hospital, date_hospitalisation) %>%   # arrange rows by hospital and then by date
     
  group_by(hospital) %>%              # group by hospital 
     
  mutate(                             # rolling average  
    mean_7day_hosp = slide_index_dbl(
      .x = new_cases,                 # the count of cases per hospital-day
      .i = date_hospitalisation,      # index on date of admission
      .f = mean,                      # use mean()                  
      .before = days(6)               # use the day and the 6 days prior
      )
  )

Đây là bộ dữ liệu mới:

Bây giờ chúng ta có thể vẽ các đường trung bình động, hiển thị dữ liệu theo nhóm bằng cách chỉ định ~ hospital tới facet_wrap() trong ggplot(). Để giải trí, chúng tôi vẽ hai biểu đồ - một biểu đồ cột geom_col() thể hiện số lượng ca mắc hàng ngày và một biểu đồ đường geom_line() thể hiện đường trung bình động của 7 ngày.

ggplot(data = grouped_roll)+
  geom_col(                       # plot daly case counts as grey bars
    mapping = aes(
      x = date_hospitalisation,
      y = new_cases),
    fill = "grey",
    width = 1)+
  geom_line(                      # plot rolling average as line colored by hospital
    mapping = aes(
      x = date_hospitalisation,
      y = mean_7day_hosp,
      color = hospital),
    size = 1)+
  facet_wrap(~hospital, ncol = 2)+ # create mini-plots per hospital
  theme_classic()+                 # simplify background  
  theme(legend.position = "none")+ # remove legend
  labs(                            # add plot labels
    title = "7-day rolling average of daily case incidence",
    x = "Date of admission",
    y = "Case incidence")

NGUY HIỂM: Nếu bạn gặp lỗi cho biết “slide() was deprecated in tsibble 0.9.0 and is now defunct. Please use slider::slide() instead.”, điều đó có nghĩa là hàm slide() từ package tsibble đang đè lên hàm slide() từ package slider. Khắc phục lỗi này bằng cách cụ thể tên package trong lệnh, ví dụ slider::slide_dbl().