---
title: "Introduction to healthdb"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Introduction to healthdb}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

------------------------------------------------------------------------

## What it does {#sec-what-it-does}

-   This package is designed for identifying disease cases from admin data for epidemiological studies. The implementation focused on code readability and re-usability. Three types of functions are included:

-   *Interactive functions* (e.g., `identify_row()`, `exclude()`, `fetch_var()`) based on filter and joins from dplyr with tweaks that fix SQL translation or add features that are not natively support by SQL. They also work for local data.frame, and some use 'data.table' package (`vignette("datatable-intro", package = "data.table")`) to speed up processing time for large data. These functions are not as flexible as `dplyr::filter()`, but they are general enough to be useful even outside health research.

-   *Call-building functions* (e.g., `build_def()`, `execute_def()`) that facilitate batch execution and re-use of case definitions. In essence, `build_def` creates codes of definitions (which is chain of the interactive functions, e.g., `define_case()`) that are not immediately ran. `execute_def` runs built definitions with different input data.

-   *Miscellaneous functions* such as computing age `compute_duration()`, collapsing records within a time range into one episode `collapse_episode()`, and more (on-going effort). Most of these functions have built-in checks signalling when things might go wrong, e.g., missing values in calculated ages.

### Motivation

In health research and surveillance, identifying diseases or events from administrative databases is often the initial step. However, crafting case-finding algorithms is a complex task. Existing algorithms, often written in SAS by experienced analysts, can be complex and difficult to decipher for the growing number of analysts trained primarily in R.

These algorithms may also affect performance if they depend on Data Step in SAS, due to a lack of translation between Data Step and SQL. This can result in SAS downloading data from a remote database to a local machine, leading to poor performance when handling large, population-based databases.

The 'healthdb' R package was created to address these challenges. It minimizes the need to download data and offers an easy-to-use interface for working with healthcare databases. It also includes capabilities not supported by 'SQL', such as matching strings by 'stringr' style regular expressions, and can compute comorbidity scores (`compute_comorbidity()`) directly on a database server. This vignette will present an example of common use cases.

------------------------------------------------------------------------

## Installation {#sec-installation}

Simply run:

```{r, eval = FALSE}
install.packages("healthdb")
```

We will need the following packages for this demo.

```{r setup, message=FALSE}
library(dplyr)
library(dbplyr)
library(lubridate)
library(glue)
library(purrr)
library(healthdb)
```

## Intended use case {#sec-intended-use-case}

Consider the case definition of substance use disorder (SUD) from [British Columbia Centre for Disease Control's Chronic Disease Dashboard](http://www.bccdc.ca/resource-gallery/Documents/Chronic-Disease-Dashboard/substance-use-disorder.pdf),

> One or more hospitalization with a substance use disorder diagnostic code, OR Two or more physician visits with a substance use disorder diagnostic code within one year.

We are going to implement this definition. First, let's make a demo data sets for the two sources:

1.  Physician claims with multiple columns of [ICD-9](https://www2.gov.bc.ca/gov/content/health/practitioner-professional-resources/msp/physicians/diagnostic-code-descriptions-icd-9) diagnostic codes

    ```{r}
# make_test_dat() makes either a toy data.frame or database table in memory with known number of rows that satisfy the query we will show later
claim_db <- make_test_dat(vals_kept = c("303", "304", "305", "291", "292", glue("30{30:59}"), glue("29{10:29}"), noise_val = c("999", "111")), type = "database")

# this is a database table
# note that in-memory SQLite database stores dates as numbers
claim_db %>% head()
    ```

2.  Hospitalization with [ICD-10](https://en.wikipedia.org/wiki/ICD-10) codes

    ```{r}
hosp_df <- make_test_dat(vals_kept = c(glue("F{10:19}"), glue("F{100:199}"), noise_val = "999"), type = "data.frame")

# this is a local data.frame/tibble
hosp_df %>% head()

# convert Date to numeric to be consistent with claim_db
hosp_df <- hosp_df %>%
  mutate(dates = julian(dates))
    ```

## Interactive functions

Let's focus on the physician claims. Extracting clients with at least two records within a year is not difficult, and involves only a few steps. The codes could look like the following using dplyr, however, it does not work because: 1. SQL does not support multiple patterns in one LIKE operation, 2. dbply currently have issue with translating n_distinct.

```{r, eval=FALSE}
## not run
claim_db %>%
  # identify the target codes
  filter(if_any(starts_with("diagx"), ~ str_like(., c("291%", "292%", "303%", "304%", "305%")))) %>%
  # each clnt has at least 2 records on different dates
  group_by(clnt_id) %>%
  # the n_distinct step is mainly for reducing computation in the next step
  filter(n_distinct(dates) >= 2) %>%
  # any two dates within one year?
  filter((max(dates) - min(dates)) <= 365)
## end
```

Here's how you could use `healthdb` to achieve these steps:

1.  Identify rows contains the target codes. Use `?identify_row` to see a list of supported matching types.

    ```{r}
result1 <- claim_db %>%
  identify_row(
    vars = starts_with("diagx"),
    match = "start",
    vals = c(291:292, 303:305)
  )
    ```

2.  Bonus: remove clients with exclusion codes

    This step is not in the substance use disorder definition, but other disease definitions often require exclusion of some ICDs that contradicts the ones of interest. Let's say we want to remove clients with code "111" here.

    We first identify "111" from the source, then exclude clients in the output from the previous step's result. `exclude()` take either a data set (via the excl argument) or expression (condition argument) as input. For the former, it performs an anti join matching on the by argument (see `dplyr::join_by()`). For the latter, it is the opposite of filter, i.e., `filter(!(some_expression))`.

    ```{r}
result2 <- result1 %>%
  exclude(
    excl = identify_row(claim_db, starts_with("diagx"), "in", "111"),
    by = "clnt_id"
  )
    ```

3.  Restrict the number of records per client

    ```{r}
result3 <- result2 %>% restrict_n(
  clnt_id = clnt_id,
  n_per_clnt = 2,
  count_by = dates,
  # here we use filter mode to remove records that failed the restriction
  mode = "filter"
)
    ```

4.  Restrict the temporal pattern of diagnoses

    `restrict_date()` supports more complicated patterns like having n diagnoses at least i days apart within j years. Note that when SQL interpret order of dates, the result could be not deterministic if there were duplicate dates within client. Therefore, a unique row id (uid) has to be supplied to get consistent result.

    ```{r}
result4 <- result3 %>% restrict_date(
  clnt_id = clnt_id,
  date_var = dates,
  n = 2,
  within = 365,
  uid = uid,
  # here we use flag mode to flag records that met the restriction instead of removing those
  mode = "flag"
)
    ```

5.  Fetch variables from other tables by matching common keys

    Up to this point, the result is only a query and have not been downloaded. Hopefully, the data has been shrunken to a manageable size for collection.

    ```{r}
# Class of result4
class(result4)

# execute query and download the result
result_df <- result4 %>% collect()

# Number of rows in source
nrow(claim_db %>% collect())

# Number of rows in the current result
nrow(result_df)
    ```

    Our data now only contains diagnoses which are probably not enough for further analyses. Let's say we want to gather client demographics such as age and sex from other sources. This certainly can be done with multiple `dplyr::left_join()` calls. Here we provide the `fetch_var()` function to make the codes more concise. Note that the input must be a named object and not from a pipe (i.e., don't do this `data %>% some_action %>% fetch_var()`). 

    ```{r}
# make two look up tables
age_tab <- data.frame(
  clnt_id = 1:50,
  age = sample(1:90, 50),
  sex = sample(c("F", "M"), 50, replace = TRUE)
)
address_tab <- data.frame(
  clnt_id = rep(1:50, 5), year = rep(2016:2020, each = 50),
  area_code = sample(0:200, 50, replace = TRUE)
)

# get year from dates for matching

result_df <- result_df %>% mutate(year = lubridate::year(as.Date(dates, origin = "1970-01-01")))

# note that keys must be present in all tables
fetch_var(result_df,
  keys = c(clnt_id, year),
  linkage = list(
    # the formula means from_table ~ get_variable
    # |clnt_id means matching on clnt_id only
    age_tab ~ c(age, sex) | clnt_id,
    address_tab ~ area_code
  )
) %>%
  select(uid, clnt_id, dates, age, sex, area_code) %>%
  head()
    ```

## Call-building functions

To complete the definition, we need to repeat the process shown above with hospitalization data. Some studies may use more than a handful of data sources to define their sample. We packed steps 1-4 in one function `define_case()`, and provide tools to perform batch execution with different data and parameters to meet those needs.

```{r}
# build the full definition of SUD
sud_def <- build_def(
  # name of definition
  def_lab = "SUD",
  # place holder names for sources
  src_labs = c("claim", "hosp"),
  def_fn = define_case, # you could alter it and supply your own function
  # below are argumets of define_case
  fn_args = list(
    # if length = 1, the single element will be use for every source
    vars = list(starts_with("diagx")),
    match = "start", # match ICD starts with vals
    vals = list(c(291:292, 303:305), glue("F{10:19}")),
    clnt_id = clnt_id,
    n_per_clnt = c(2, 1),
    date_var = dates,
    within = c(365, NULL),
    uid = uid,
    mode = "flag"
  )
)

sud_def
```

Let's look inside the fn_call list column. Two calls of `define_case()` have been made with different parameters. The data arguments are left empty on purpose for re-usability. For example, you may want to repeat the analysis with data from different regions or study periods.

```{r}
sud_def$fn_call
```

Executing the definition is simply a call of `execute_def()`. If verbose option is not turned off by `options(healthdb.verbose = FALSE)`, the output message will explain what has been done. You could append multiple `build_def()` outputs together and execute them all at once. Definition and source labels will be added to the result to identify outputs from different calls.

```{r}
# execute the definition
result_list <- sud_def %>%
  execute_def(with_data = list(
    claim = claim_db,
    hosp = hosp_df
  ))
```

Let's check the results!

```{r}
# view the results
purrr::walk(result_list, ~ head(.) %>% print())
```

At this point, the result from the claim database (`result[[1]]`) has not been collected locally. You could collect it manually, do further filtering, and then combine with the result from hospitalization data in any way you want. If you just need a simple row bind, we have `bind_source()` with convenient naming feature.

```{r}
bind_source(result_list,
  # output_name = c(names in the list elements)
  src = "src",
  uid = "uid",
  clnt_id = "clnt_id",
  flag = c("flag_restrict_date", NA),
  # force_proceed is needed to collect remote tables to local memory
  force_proceed = TRUE
)
```

`pool_case()` goes a few steps further than row bind. It also filters records with valid flags and can summarize by client/group. Since we had to decide which variables to be summarized in advance, the output may not be flexible enough to meet your needs.

```{r}
pool_case(result_list,
  def = sud_def,
  # your could skip summary with output_lvl = "raw"
  output_lvl = "clnt",
  # include records only from sources having valid records, see function documentation for more detail and other options
  include_src = "has_valid",
  force_proceed = TRUE
)
```
