Follow-up courses

We offer the following follow-up courses within the theme Core Methods of the Public Health Sciences Course Program:

Timetable

Today

Course structure and material

• Lectures: HTML slides
• Exercises: .R and .qmd files
• Assessment: Published report on your GitHub account
• GitHub repository: https://github.com/ISPMBern/basic-statistics-and-projects-in-R
• GitHub Pages site: https://ispmbern.github.io/basic-statistics-and-projects-in-R
• Evaluation: https://www.tiny.cc/phs23 (or QR code), TAN: SPQPL
• Interactive!

Installing R, RStudio, and Git

Mandatory steps (should be completed by now)

1. Download and install R.
2. Download and install RStudio.
3. Install Git using the following instructions.
4. If you don’t already have one, don’t forget to create a GitHub account.

Optional steps (can be done today)

1. Make sure RStudio knows about Git by following the corresponding section here.
2. Install the usethis package for R using the following command: install.packages("usethis")
3. Set up Git using the following command: usethis::use_git_config(user.name = "Jane Doe", user.email = "jane@example.org")
4. Generate a personal access token (PAT) and store your PAT as described here.

What is R?

• R is a programming language for statistical computing and graphics that first appeared in 1993.
• R is the open source implementation of the S language, which was developed by Bell laboratories in the 1970’s.
• The aim of the S language, as expressed by John Chambers, is “to turn ideas into software, quickly and faithfully”.
• The statisticians Ross Ihaka and Robert Gentleman developed R at the University of Auckland, New Zealand.
• In 1995, the statistician Martin Mächler from ETH Zurich convinced Ihaka and Gentleman to make R free and open source under the GNU General Public License.
• R is both open source and open development.

Why R?

• Free and open source (popular in academic research)
• High-level programming language designed for statistical computing (popular in computational biology, bioinformatics, epidemiology and public health sciences)
• Reproducibility
• Powerful and flexible - especially for data wrangling and visualization
• Extensive add-on software (packages)
• Strong community

(Why not R?)

• Little centralized support, relies on online community and package developers
• Can be cumbersome to update
• Slower than more traditional programming languages (C/C++, Python)

R vs. RStudio

R is a programming language that runs computations. RStudio is an integrated development environment (IDE) that provides an interface by adding many convenient features and tools, e.g.:

• autocomplete, syntax, and spell checking functionality,
• easier working with file paths (via projects and the here package),
• integration with version control systems (e.g., GitHub, SVN).

Note that RStudio (the company) recently changed it’s name to Posit. RStudio (the IDE) remains unchanged.

RStudio cheatsheet

Use the cheatsheet to find your way in RStudio: https://github.com/rstudio/cheatsheets/blob/main/rstudio-ide.pdf

Object classes

All the objects stored in R have a class which tells R how to handle the object. There are many possible classes, but common ones include:

Other object classes are array, data.frame, list, and tibble (similar to data.frame).

Writing your own functions

Every function in R has three basic parts: a name, a body of code, and a set of arguments. To make your own function, you need to replicate these parts and store them in an R object, which you can do with the function function().

The R ecosystem

The build-in functionality of Base R can be expaned with packages that others have developed and published.

The Comprehensive R Archive Network (CRAN) has been the main source of R packages. Nowadays, GitHub also contains many packages and has arguably become the primary location for package development. Some packages, e.g., tidyverse, are so-called meta-packages - they load a collection of other packages.

Install new packages as follows:

install.packages("packagename")

Packages must be loaded each R session to give access to their functionality:

library(packagename)

Further resources

There are numerous online tutorials and books on R, RStudio with specific applications to epidemiology, public health, and data science:

• The R Manuals
• Hands-On Programming with R
• The Epidemiologist R Handbook
• Introduction to R for Public Health Researchers
• R for Data Science

Why use projects in R(Studio)?

The use of projects (.Rproj files) is fundamental to organized coding and project management. There are four main reasons why you should use projects essentially 100% of the time while using RStudio:

• It can take less than 30 seconds to set up.
• It keeps all relevant files in the same place.
• It sets the working directory, so you can use relative paths.
• It allows for version control (see tomorrow).

Folder structures

Using consistent folder structures across projects helps you work more efficiently. An example of the folder structure for a project looks like this:

• R
  • 00_main.R
  • 01_cleaning.R
  • 02_analysis.R
  • 03_plotting.R
• data
  • processed
    • cleaned_data.csv
    • processed_data.rds
  • raw
    • original_data.csv
    • spreadsheet_data.xlsx
• output
  • figures
    • 01_figure.png
    • 02_figure.pdf
  • tables
    • 01_table.csv
    • 02_table.rds
• products
  • manuscript
    • manuscript.docx
    • manuscript.html
    • manuscript.pdf
    • manuscript.qmd
  • report
    • report.html
    • report.qmd
  • slides
    • slides.html
    • slides.qmd
• .gitignore
• README.md
• project-template.Rproj

Naming files

Bad examples

• myabstract.docx
• Long file name using “spaces” & punctuation!.xlsx
• figures 1.png
• SDEFI7_jknsfol.txt

Better

• 2023-02-15_abstract-conference-X.docx
• still-long-but-no-punctuation-or-spaces.xlsx
• fig01_scatter-mpg-vs-vol.png
• more-meaningful-name.txt

Naming R files

Analysis scripts can get looooooooong… Don’t be afraid to break them up into smaller chunks.

Use sequential numbers and descriptive names, e.g.:

• 01_cleaning.R cleans your data,
• 02_analysis.R performs your analysis,
• 03_plotting.R plots the results of your analysis.

Sequential numbers allow you to sort the files according to the sequence in which you run them.

Descriptive names inform you of what is actually in there.

You can use a main or master file (00_main.R) to run all other files and create a reproducible analysis.

Also see How to name files from Jennifer Bryan and The tidyverse style guide.

Required packages for the next exercises

To work on the next exercises, you have to install the following packages:

• usethis - Workflow package
• gitcreds - Queries Git credentials from R
• here - Easy file referencing
• tidyverse - A set of packages
• medicaldata - Medical data sets
• cowplot - Features to create publication-quality figures

Simply type install.packages("packagename"), but RStudio will ask you about it as well if you want to load a package that you haven’t installed yet.

What is the tidyverse?

What is the tidyverse?

Getting data into R

R has a wide range of tools for importing data.

Getting data into R

R has a wide range of tools for importing data.

Getting data into R

R has a wide range of tools for importing data.

Data sometimes already exists in R…

Many published datasets already exist in R, either in the basic installation or via a package

data(mtcars)
data(iris) # very popular in examples

From packages, e.g. medicaldata¹

install.packages("medicaldata")
library(medicaldata)
strep_tb

# A tibble: 107 × 13
   patient_id arm     dose_strep_g dose_PAS_g gender baseline_condition
   <chr>      <fct>          <dbl>      <dbl> <fct>  <fct>             
 1 0001       Control            0          0 M      1_Good            
 2 0002       Control            0          0 F      1_Good            
 3 0003       Control            0          0 F      1_Good            
 4 0004       Control            0          0 M      1_Good            
 5 0005       Control            0          0 F      1_Good            
 6 0006       Control            0          0 M      1_Good            
 7 0007       Control            0          0 F      1_Good            
 8 0008       Control            0          0 M      1_Good            
 9 0009       Control            0          0 F      2_Fair            
10 0010       Control            0          0 M      2_Fair            
# ℹ 97 more rows
# ℹ 7 more variables: baseline_temp <fct>, baseline_esr <fct>,
#   baseline_cavitation <fct>, strep_resistance <fct>, radiologic_6m <fct>,
#   rad_num <dbl>, improved <lgl>

Getting data into R in practice

We place the dataset in the appropriate folder (01_original_data) and read it in with the appropriate function (e.g. read_csv)

library(readr)
data <- read_csv(here("data", "raw", "MyData.csv"))

Depending on the file, you might need read_csv2, which is configured for e.g. German environments where CSVs are actually semi-colon (;) separated, because the comma is used in numbers . . .

In base R

data <- read.csv(here("data", "raw", "MyData.csv"))

Virtually identical… readr is slightly faster and automatically converts some variable types

Getting data into R in practice

More common(?): excel files…

library(readxl) # informal tidyverse member
data <- read_xlsx(here("data", "raw", "MyData.xlsx"))

str(strep_tb)

tibble [107 × 13] (S3: tbl_df/tbl/data.frame)
 $ patient_id         : chr [1:107] "0001" "0002" "0003" "0004" ...
 $ arm                : Factor w/ 2 levels "Streptomycin",..: 2 2 2 2 2 2 2 2 2 2 ...
 $ dose_strep_g       : num [1:107] 0 0 0 0 0 0 0 0 0 0 ...
 $ dose_PAS_g         : num [1:107] 0 0 0 0 0 0 0 0 0 0 ...
 $ gender             : Factor w/ 2 levels "F","M": 2 1 1 2 1 2 1 2 1 2 ...
 $ baseline_condition : Factor w/ 3 levels "1_Good","2_Fair",..: 1 1 1 1 1 1 1 1 2 2 ...
 $ baseline_temp      : Factor w/ 4 levels "1_98-98.9F","2_99-99.9F",..: 1 3 1 1 2 3 2 2 2 4 ...
 $ baseline_esr       : Factor w/ 4 levels "1_0-10","2_11-20",..: 2 2 3 3 3 3 3 3 3 3 ...
 $ baseline_cavitation: Factor w/ 2 levels "no","yes": 2 1 1 1 1 1 2 2 2 2 ...
 $ strep_resistance   : Factor w/ 3 levels "1_sens_0-8","2_mod_8-99",..: 1 1 1 1 1 1 1 1 1 1 ...
 $ radiologic_6m      : Factor w/ 6 levels "6_Considerable_improvement",..: 1 2 2 2 2 1 2 2 2 2 ...
 $ rad_num            : num [1:107] 6 5 5 5 5 6 5 5 5 5 ...
 $ improved           : logi [1:107] TRUE TRUE TRUE TRUE TRUE TRUE ...

XLSX vs CSV vs database (e.g. REDCap)

Watch Darren Dahly open an Excel file

What do you think this is?

[1] "\xfc" "\xe4" "\xe9" "\xe0"

Dealing with special characters

[1] "\xfc" "\xe4" "\xe9" "\xe0"

[1] "ü" "ä" "é" "à"

# base
read.csv("path/to/file.csv", fileEncoding = "UTF-8")
# tidyverse (readr)
read_csv("path/to/file.csv", fileEncoding = "UTF-8",
         locale = locale(encoding = "UTF-8"))

Solution

library(readr)
dat <- read_csv("data/raw/insurance_with_date.csv")
str(dat)

spc_tbl_ [1,338 × 9] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
 $ X       : num [1:1338] 1 2 3 4 5 6 7 8 9 10 ...
 $ age     : num [1:1338] 59 24 28 22 60 38 51 44 47 29 ...
 $ sex     : chr [1:1338] "male" "female" "female" "male" ...
 $ bmi     : num [1:1338] 31.8 22.6 25.9 25.2 36 ...
 $ children: num [1:1338] 2 0 1 0 0 3 0 0 1 2 ...
 $ smoker  : chr [1:1338] "no" "no" "no" "no" ...
 $ region  : chr [1:1338] "southeast" "southwest" "northwest" "northwest" ...
 $ charges : num [1:1338] 13086 2574 4411 2321 13435 ...
 $ date    : Date[1:1338], format: "2001-01-15" "2001-01-17" ...
 - attr(*, "spec")=
  .. cols(
  ..   X = col_double(),
  ..   age = col_double(),
  ..   sex = col_character(),
  ..   bmi = col_double(),
  ..   children = col_double(),
  ..   smoker = col_character(),
  ..   region = col_character(),
  ..   charges = col_double(),
  ..   date = col_date(format = "")
  .. )
 - attr(*, "problems")=<externalptr> 

dat2 <- read.csv("data/raw/insurance_with_date.csv")
str(dat2)

'data.frame':   1338 obs. of  9 variables:
 $ X       : int  1 2 3 4 5 6 7 8 9 10 ...
 $ age     : int  59 24 28 22 60 38 51 44 47 29 ...
 $ sex     : chr  "male" "female" "female" "male" ...
 $ bmi     : num  31.8 22.6 25.9 25.2 36 ...
 $ children: int  2 0 1 0 0 3 0 0 1 2 ...
 $ smoker  : chr  "no" "no" "no" "no" ...
 $ region  : chr  "southeast" "southwest" "northwest" "northwest" ...
 $ charges : num  13086 2574 4411 2321 13435 ...
 $ date    : chr  "2001-01-15" "2001-01-17" "2001-01-22" "2001-01-29" ...

Keeping information together - lists

For when the objects have different lengths…

lst <- list(letter = "a",
     numbers = rnorm(10),
     data = df)
lst

$letter
[1] "a"

$numbers
 [1]  1.6767047 -1.8662187  1.0433395 -0.6191508  0.4764069  2.3812988
 [7] -0.7050137  1.1650151  0.6387094 -0.6682568

$data
     sex    fct num   lgl
1   male   male   1  TRUE
2 female female   2 FALSE
3 female female   3 FALSE
4   male   male   4  TRUE

Getting elements out again is the same for both

df$sex

[1] "male"   "female" "female" "male"  

lst$letter

[1] "a"

Modifying data

mutate is your friend

strep_tb |> 
  mutate(dose_strep_g = dose_strep_g + 2)

Very simple in base too…

strep_tb$dose_strep_g <- strep_tb$dose_strep_g + 2

strep_tb |> 
  mutate(dose_strep_g = dose_strep_g + 2, 
         control = arm == "Control")

within(strep_tb, {
  dose_strep_g <- dose_strep_g + 2
  control <- arm == "Control"
})

# A tibble: 107 × 14
   patient_id arm     dose_strep_g dose_PAS_g gender baseline_condition
   <chr>      <fct>          <dbl>      <dbl> <fct>  <fct>             
 1 0001       Control            2          0 M      1_Good            
 2 0002       Control            2          0 F      1_Good            
 3 0003       Control            2          0 F      1_Good            
 4 0004       Control            2          0 M      1_Good            
 5 0005       Control            2          0 F      1_Good            
 6 0006       Control            2          0 M      1_Good            
 7 0007       Control            2          0 F      1_Good            
 8 0008       Control            2          0 M      1_Good            
 9 0009       Control            2          0 F      2_Fair            
10 0010       Control            2          0 M      2_Fair            
# ℹ 97 more rows
# ℹ 8 more variables: baseline_temp <fct>, baseline_esr <fct>,
#   baseline_cavitation <fct>, strep_resistance <fct>, radiologic_6m <fct>,
#   rad_num <dbl>, improved <lgl>, control <lgl>

Working with strings: stringr

The stringr package contains functions specifically for working with strings.

Most functions start with str_.

library(stringr)
txt <- "A siLLy exAmple "

str_to_lower(txt)

[1] "a silly example "

str_to_sentence(txt)

[1] "A silly example "

str_to_title(txt)

[1] "A Silly Example "

str_to_upper(txt)

[1] "A SILLY EXAMPLE "

str_length(txt)

[1] 16

str_count(txt, "\\w+")

[1] 3

Working with strings: stringr

Remove white space

str_squish(txt)

[1] "A siLLy exAmple"

str_sub(txt, 3)

[1] "siLLy exAmple "

str_sub(txt, 3, -4)

[1] "siLLy exAmp"

word(txt, 2)

[1] "siLLy"

str_replace(txt, "e", "XX")

[1] "A siLLy XXxAmple "

str_replace_all(txt, "e", "XX")

[1] "A siLLy XXxAmplXX "

Working with strings: stringr

Detect a substring

str_detect(txt, "x")

[1] TRUE

str_detect(txt, "z")

[1] FALSE

str_split(txt, " ")

[[1]]
[1] "A"       "siLLy"   "exAmple" ""       

strep_tb |> 
  mutate(txt = as.character(baseline_condition),
         upper = str_to_upper(txt),
         good = str_detect(txt, "Good"),
         no_number = str_replace(txt, "^[[:digit:]]_", "")) |> 
  select(txt:no_number) |> unique()

# A tibble: 3 × 4
  txt    upper  good  no_number
  <chr>  <chr>  <lgl> <chr>    
1 1_Good 1_GOOD TRUE  Good     
2 2_Fair 2_FAIR FALSE Fair     
3 3_Poor 3_POOR FALSE Poor     

Working with factors: forcats

The forcats package contains functions specifically for working with factors.

Functions (almost) all begin with fct_.

fac <- strep_tb$baseline_condition[c(1, 15, 29)]
fac

[1] 1_Good 2_Fair 3_Poor
Levels: 1_Good 2_Fair 3_Poor

library(forcats)
fct_rev(fac)

[1] 1_Good 2_Fair 3_Poor
Levels: 3_Poor 2_Fair 1_Good

# factor(fac, 
#        levels = rev(levels(fac)))

Working with factors: forcats

Change level names

fct_recode(fac,
           # new = "old"
           Good = "1_Good",
           Fair = "2_Fair",
           Poor = "3_Poor")

[1] Good Fair Poor
Levels: Good Fair Poor

# factor(as.character(fac),
#        levels = levels(fac),
#        labels = c("Good", "Fair", "Poor"))

fct_relabel(fac,
            str_replace, # function to use
            # additional arguments to the function
            pattern = "^[123]_", # regex for "starts with 1, 2 or 3 and is followed by _"
            replacement = "" # replace with nothing
            ) 

[1] Good Fair Poor
Levels: Good Fair Poor

# factor(as.character(fac),
#        levels = levels(fac),
#        labels = gsub("^[123]_", "", levels(fac)))

Working with factors: forcats

Within mutate

strep_tb |> 
  mutate(baseline_condition_new = fct_relabel(baseline_condition,
                                              str_replace, 
                                              pattern = "^[[:digit:]]_", # "any leading digit followed by underscore"
                                              replacement = "")) |> 
  select(baseline_condition, baseline_condition_new) |> 
  str()

tibble [107 × 2] (S3: tbl_df/tbl/data.frame)
 $ baseline_condition    : Factor w/ 3 levels "1_Good","2_Fair",..: 1 1 1 1 1 1 1 1 2 2 ...
 $ baseline_condition_new: Factor w/ 3 levels "Good","Fair",..: 1 1 1 1 1 1 1 1 2 2 ...

strep_tb |> 
  mutate(across(where(is.factor), 
                ~ fct_relabel(.x,
                              str_replace, 
                              pattern = "^[[:digit:]]_", 
                              replacement = ""))) |> 
  select(where(is.factor)) |> 
  str()

tibble [107 × 8] (S3: tbl_df/tbl/data.frame)
 $ arm                : Factor w/ 2 levels "Streptomycin",..: 2 2 2 2 2 2 2 2 2 2 ...
 $ gender             : Factor w/ 2 levels "F","M": 2 1 1 2 1 2 1 2 1 2 ...
 $ baseline_condition : Factor w/ 3 levels "Good","Fair",..: 1 1 1 1 1 1 1 1 2 2 ...
 $ baseline_temp      : Factor w/ 4 levels "98-98.9F","99-99.9F",..: 1 3 1 1 2 3 2 2 2 4 ...
 $ baseline_esr       : Factor w/ 4 levels "0-10","11-20",..: 2 2 3 3 3 3 3 3 3 3 ...
 $ baseline_cavitation: Factor w/ 2 levels "no","yes": 2 1 1 1 1 1 2 2 2 2 ...
 $ strep_resistance   : Factor w/ 3 levels "sens_0-8","mod_8-99",..: 1 1 1 1 1 1 1 1 1 1 ...
 $ radiologic_6m      : Factor w/ 6 levels "Considerable_improvement",..: 1 2 2 2 2 1 2 2 2 2 ...

Working with dates: lubridate

lubridate provides a comprehensive set of functions for working with dates and date-times

Dates come in many formats, lubridate handles them easily

library(lubridate)
ymd("2023-01-20")

[1] "2023-01-20"

ymd_hm("2023-01-20 10:15")

[1] "2023-01-20 10:15:00 UTC"

dmy("20 January 2023")

[1] "2023-01-20"

mdy("January 20 2023")

[1] "2023-01-20"

mdy("Januar 20 2023")

[1] NA

mdy("1 20 2023")

[1] "2023-01-20"

mdy("1 20 2023")

[1] "2023-01-20"

dmy("the 1st of May 2023 was a Monday")

[1] "2023-05-01"

Working with dates: base

In base R, it’s not so easy… (see ?strptime for details)

as.Date("2023-01-20")

[1] "2023-01-20"

as.Date("20 Jan 2023", format = "%d %B %Y")

[1] "2023-01-20"

as.Date("01 20 2023", format = "%m %d %Y")

[1] "2023-01-20"

as.POSIXct("2023-01-20 10:15")

[1] "2023-01-20 10:15:00 CET"

Very specific to system settings (language)…

Unless you have e.g. “2023-01-20 10:15”, stick with lubridate…

Merging datasets (joins)

(a <- data.frame(id = 1:3, 
                 v1 = letters[1:3]))

  id v1
1  1  a
2  2  b
3  3  c

(b <- data.frame(id = 2:4, 
                 v2 = LETTERS[1:3]))

  id v2
1  2  A
2  3  B
3  4  C

a |> left_join(b)

  id v1   v2
1  1  a <NA>
2  2  b    A
3  3  c    B

a |> right_join(b)

  id   v1 v2
1  2    b  A
2  3    c  B
3  4 <NA>  C

a |> inner_join(b)

  id v1 v2
1  2  b  A
2  3  c  B

Merging datasets (joins)

a |> full_join(b)

  id   v1   v2
1  1    a <NA>
2  2    b    A
3  3    c    B
4  4 <NA>    C

a |> full_join(b, 
               by = join_by(id_a == id_b) # new syntax!
               # by = c("id_a" = "id_b") # older syntax
               )

merge(a, b, all.x = TRUE,  all.y = FALSE, by.x, by.y)  # left
merge(a, b, all.x = FALSE, all.y = TRUE,  by.x, by.y)  # right
merge(a, b, all.x = FALSE, all.y = FALSE, by.x, by.y)  # inner
merge(a, b, all.x = TRUE,  all.y = TRUE,  by.x, by.y)  # full

Summarizing data

At some point, you will have to create summary data. dplyr can help with that too. summarize is the appropriate function.

strep_tb |> 
  summarize(n = n(),
            min = min(rad_num),
            median = median(rad_num),
            mean = mean(rad_num),
            max = min(rad_num),
            )

# A tibble: 1 × 5
      n   min median  mean   max
  <int> <dbl>  <dbl> <dbl> <dbl>
1   107     1      5  3.93     1

strep_tb |> 
  summarize(n = n(),
            across(c(rad_num, dose_strep_g), 
                   list(min = ~ min(.x, na.rm = TRUE), 
                        mean = mean, 
                        median = median, 
                        max = max), 
                   .names = "{.col}_{.fn}"),
            )

# A tibble: 1 × 9
      n rad_num_min rad_num_mean rad_num_median rad_num_max dose_strep_g_min
  <int>       <dbl>        <dbl>          <dbl>       <dbl>            <dbl>
1   107           1         3.93              5           6                0
# ℹ 3 more variables: dose_strep_g_mean <dbl>, dose_strep_g_median <dbl>,
#   dose_strep_g_max <dbl>

Summarizing data

What about grouped summaries? Two options… group_by()…

strep_tb |> 
  group_by(arm) |> 
  summarize(n = n(),
            min = min(rad_num),
            median = median(rad_num),
            mean = mean(rad_num),
            max = min(rad_num),
            )

… or .by (new syntax)

strep_tb |> 
  summarize(n = n(),
            min = min(rad_num),
            median = median(rad_num),
            mean = mean(rad_num),
            max = min(rad_num),
            .by = arm
            )

# A tibble: 2 × 6
  arm              n   min median  mean   max
  <fct>        <int> <dbl>  <dbl> <dbl> <dbl>
1 Control         52     1      3  3.13     1
2 Streptomycin    55     1      6  4.67     1

Publication type tables

gtsummary is my package of choice

library(gtsummary)
strep_tb_lab |> 
  select(arm, dose_strep_g, baseline_temp,
         rad_num, improved) |> 
  tbl_summary(by = arm, 
              type = c(
                rad_num ~ "continuous"
              )) |> 
  add_overall()