monthlytemporalVariability.Rmd

---
title: 'Monthly Temporal Variability in Pediatric population: SITEID'
output:
  html_document:
    code_folding: show
    toc: yes
    toc_float: yes
  pdf_document:
    toc: yes
  word_document:
    toc: yes
---

```{r libraries, message=FALSE, warning=FALSE, results='hide'}
# Install packages
paket <- function(pak){
  new_pak <- pak[!(pak %in% rownames(installed.packages()))]
  if (length(new_pak)) 
    install.packages(new_pak, dependencies = TRUE,repos="https://cloud.r-project.org/")
  sapply(pak, library, character.only = TRUE)
}

listOfPackages <- c("EHRtemporalVariability", "dplyr", "tidyr", "tidyverse", "tableone", "plotly", "DT")
paket(listOfPackages)
```

```{r}
theme_set(theme_bw() + 
            theme(legend.title = element_blank(),
                  panel.grid.minor = element_blank()))
```


## Data read-in

Read the functions and data needed for the analysis. 
```{r suppData, message=FALSE}
source("R/readInputFiles.R")
source("R/summaryTableOne.R")
source("R/tempVariability.R")
source("R/countsPerPatientMonth.R")
source("R/countsPerAgeGroup.R")
mentalDisorder <- read.delim('public-data/icd10_categories_mentalDisorders.txt')
icdCategory <- read.delim("public-data/icd10Codes.txt", header = FALSE, colClasses = "character")
```

Read the input files.

```{r message=FALSE, warning=FALSE}
### Read the CSV input files
#change to the path where your phase 2.1 data files are located
myinputFiles <- "/4ceData/Input/"
files <- readInputFiles( path      = myinputFiles, 
                         separator = ",",
                         skip      = 0, 
                         verbose   = FALSE )
  
### Extract the patient summary and observation information. 
demo_raw <- files[[1]]
obs_raw <- files[[2]]
clinical_raw <- files[[3]]

### Select only pediatric subset
demo_raw     <- demo_raw[ demo_raw$age_group %in% c("00to02", "03to05", "06to11","12to17"), ]
obs_raw      <- obs_raw[ obs_raw$patient_num %in% demo_raw$patient_num, ]
clinical_raw <- clinical_raw[ clinical_raw$patient_num %in% demo_raw$patient_num, ]

### Determine the obfuscation threshold (FALSE if no obfuscation, the numeric value of the obfuscation threshold if any)
obfuscation =  FALSE
#obfuscation = 3

### Determine if you don't want to visualize counts
absoluteValue = TRUE

### Determine the analysis period
period = "month"
```

```{r}
counts_all <- countsPerPatientMonth( observationFile = obs_raw, clinicalFile = clinical_raw, type = "diagnosis", obfuscation = obfuscation)
counts_all[[1]] %>% ggplot() +
  geom_line(aes(x = year_month,y = distinct_patients, group=1), color = "red",linetype = "dotted") +
  labs(x = "calendar months", y = "number patients consulting for condition")+
  theme_bw()+
  theme(strip.text.x = element_text(size = 7), text = element_text(size=9), 
        axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))
```

# Mental Disorders
ICD10 codes extracted from: The ICD-10 Classification of Mental and Behavioural Disorders Clinical descriptions and diagnostic guidelines by the World Health Organization https://www.who.int/classifications/icd/en/bluebook.pdf

Suidice related codes added from: https://www.cdc.gov/nchs/data/nhsr/nhsr108.pdf 

```{r}
DT::datatable( mentalDisorder, rownames = FALSE )
```

```{r  message=FALSE, warning=FALSE}
probmaps_md <- tempVariability( VarSelection = "mental", 
                     observations = obs_raw,
                     clinical    = clinical_raw, 
                     demographics = demo_raw, 
                     period = period, 
                     obfuscation = obfuscation, 
                     absoluteValue = absoluteValue)

plotDataTemporalMap(
      dataTemporalMap =  probmaps_md[["mentaldescriptionShort"]],
      startValue = 1,
      endValue = 20,
      startDate = "2020-04-01",
      colorPalette    = "Magma", 
      absolute = absoluteValue)
```

# Aggregated by PheWAS category
PheCode mapping file downloaded from: https://phewascatalog.org/phecodes_icd10 


```{r  message=FALSE, warning=FALSE}
probmaps_pc <- tempVariability( VarSelection = "phecodes", 
                     observations = obs_raw,
                     clinical    = clinical_raw, 
                     demographics = demo_raw, 
                     period = period, 
                     obfuscation = obfuscation, 
                     absoluteValue = absoluteValue)

plotDataTemporalMap(
      dataTemporalMap =  probmaps_pc[["phecodesdescriptionShort"]],
      startValue = 1,
      endValue = 20,
      startDate = "2020-04-01",
      colorPalette    = "Magma", 
      absolute = absoluteValue)
```

# All ICD codes

```{r  message=FALSE, warning=FALSE}
probmaps_all <- tempVariability( VarSelection = "all", 
                     observations = obs_raw,
                     clinical    = clinical_raw, 
                     demographics = demo_raw, 
                     period = period, 
                     obfuscation = obfuscation, 
                     absoluteValue = absoluteValue)
 
plotDataTemporalMap(
       dataTemporalMap =  probmaps_all$allconcept_code,
       startValue = 1,
       endValue = 20,
       startDate = "2020-03-01",
       colorPalette    = "Magma", 
       absolute = absoluteValue)
```

# Higher general categories


```{r}
DT::datatable( icdCategory, rownames = FALSE )
```

```{r}
probmaps_category <- tempVariability( VarSelection = "category", 
                     observations = obs_raw,
                     clinical    = clinical_raw, 
                     demographics = demo_raw, 
                     period = period, 
                     obfuscation = obfuscation, 
                     absoluteValue = absoluteValue)

plotDataTemporalMap(
       dataTemporalMap =  probmaps_category$categoryDescription,
       startValue = 1,
       endValue = 20,
       startDate = "2020-03-01",
       colorPalette    = "Magma", 
       absolute = absoluteValue)
```

```{r}
counts <- countsPerAgeGroup( observationFile = obs_raw, clinicalFile = clinical_raw, demogFile = demo_raw, obfuscation = FALSE)
```

# Medication

```{r}
countsMeds <- countsPerPatientMonth( observationFile = obs_raw, clinicalFile = clinical_raw, type = "medication", obfuscation = obfuscation)
countsMeds[[1]]%>% ggplot() +
  geom_line(aes(x = year_month,y = distinct_patients, group=1), color = "red",linetype = "dotted") +
  labs(x = "calendar months", y = "number patients consulting for condition")+
  theme_bw()+
  theme(strip.text.x = element_text(size = 7), text = element_text(size=9), 
        axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))
```


```{r  message=FALSE, warning=FALSE}
probmaps_meds <- tempVariability( VarSelection = "medication", 
                     observations = obs_raw,
                     clinical    = clinical_raw, 
                     demographics = demo_raw, 
                     period = period, 
                     obfuscation = obfuscation, 
                     absoluteValue = absoluteValue)

plotDataTemporalMap(
      dataTemporalMap =  probmaps_meds[["medicationconcept_code"]],
      startValue = 1,
      endValue = 20,
      colorPalette    = "Magma", 
      absolute = absoluteValue)
```

# Save the RData file. 
If the RData file is not correctly generated, you can try to change the extension to .rda, uncomment the last line of code in this chunk.
```{r}
probMaps <- c(probmaps_md, probmaps_all, probmaps_pc, probmaps_category, probmaps_meds)
igtProjs <- sapply ( probMaps, estimateIGTProjection )
countsDf <- c( counts_all, counts, countsMeds)
save(probMaps, igtProjs, countsDf, file = paste0("myExport_", period, ".RData"))
#save(probMaps, igtProjs, file = paste0("myExport_", period, ".rda"))
```


# R session information
```{r}
sessionInfo()
```