Michael's and Christian's Blog

Category: Programming

  • DuckDB: Quacking SQL

    Lost in Translation between R and Python 8

    This is the next article in our series “Lost in Translation between R and Python”. The aim of this series is to provide high-quality R and Python 3 code to achieve some non-trivial tasks. If you are to learn R, check out the R tab below. Similarly, if you are to learn Python, the Python tab will be your friend.

    DuckDB

    DuckDB is a fantastic in-process SQL database management system written completely in C++. Check its official documentation and other blogposts like this to get a feeling of its superpowers. It is getting better and better!

    Some of the highlights:

    • Easy installation in R and Python, made possible via language bindings.
    • Multiprocessing and fast.
    • Allows to work with data bigger than RAM.
    • Can fire SQL queries on R and Pandas tables.
    • Can fire SQL queries on (multiple!) csv and/or Parquet files.
    • Quacks Apache Arrow.

    Installation

    DuckDB is super easy to install:

    • R: install.packages("duckdb")
    • Python: pip install duckdb

    Additional packages required to run the code of this post are indicated in the code.

    A first query

    Let’s start by loading a dataset, initializing DuckDB and running a simple query.

    The dataset we use here contains information on over 20,000 sold houses in Kings County. Along with the sale price, different features describe the size and location of the properties. The dataset is available on OpenML.org with ID 42092.

    library(OpenML)
library(duckdb)
library(tidyverse)

# Load data
df <- getOMLDataSet(data.id = 42092)$data

# Initialize duckdb, register df and materialize first query
con = dbConnect(duckdb())
duckdb_register(con, name = "df", df = df)
con %>% 
  dbSendQuery("SELECT * FROM df limit 5") %>% 
  dbFetch()
    import duckdb
import pandas as pd
from sklearn.datasets import fetch_openml

# Load data
df = fetch_openml(data_id=42092, as_frame=True)["frame"]

# Initialize duckdb, register df and fire first query
# If out-of-RAM: duckdb.connect("py.duckdb", config={"temp_directory": "a_directory"})
con = duckdb.connect()
con.register("df", df)
con.execute("SELECT * FROM df limit 5").fetchdf()
    Result of first query (from R)

    Average price per grade

    If you like SQL, then you can do your data preprocessing and simple analyses with DuckDB. Here, we calculate the average house price per online grade (the higher the grade, the better the house).

    query <- 
  "
  SELECT AVG(price) avg_price, grade 
  FROM df 
  GROUP BY grade
  ORDER BY grade
  "
avg <- con %>% 
  dbSendQuery(query) %>% 
  dbFetch()

avg

    # Average price per grade
query = """
  SELECT AVG(price) avg_price, grade 
  FROM df 
  GROUP BY grade
  ORDER BY grade
  """
avg = con.execute(query).fetchdf()
avg
    R output

    Highlight: queries to files

    The last query will be applied directly to files on disk. To demonstrate this fantastic feature, we first save “df” as a parquet file and “avg” as a csv file. 

    write_parquet(df, "housing.parquet")
write.csv(avg, "housing_avg.csv", row.names = FALSE)

    # Save df and avg to different file types
df.to_parquet("housing.parquet")  # pyarrow=7
avg.to_csv("housing_avg.csv", index=False)

    Let’s load some columns of “housing.parquet” data, but only rows with grades having an average price of one million USD. Agreed, that query does not make too much sense but I hope you get the idea…😃

    # "Complex" query
query2 <- "
  SELECT price, sqft_living, A.grade, avg_price
  FROM 'housing.parquet' A
  LEFT JOIN 'housing_avg.csv' B
  ON A.grade = B.grade
  WHERE B.avg_price > 1000000
  "

expensive_grades <- con %>% 
  dbSendQuery(query2) %>% 
  dbFetch()

head(expensive_grades)

# dbDisconnect(con)
    # Complex query
query2 = """
  SELECT price, sqft_living, A.grade, avg_price
  FROM 'housing.parquet' A
  LEFT JOIN 'housing_avg.csv' B
  ON A.grade = B.grade
  WHERE B.avg_price > 1000000
  """
expensive_grades = con.execute(query2).fetchdf()
expensive_grades

# con.close()
    R output

    Last words

    • DuckDB is cool!
    • If you have strong SQL skills but do not know R or Python so well, this is a great way to get used to those programming languages.
    • If you are unfamiliar to SQL but like R and/or Python, you can use DuckDB for a while and end up being an SQL addict.
    • If your analysis involves combining many large files during preprocessing, then you can try the trick shown in the last example of this post.

    The Python notebook and R code can be found at:

  • Avoid loops in R! Really?

    It must have been around the year 2000, when I wrote my first snipped of SPLUS/R code. One thing I’ve learned back then:

    Loops are slow. Replace them with

    1. vectorized calculations or
    2. if vectorization is not possible, use sapply() et al.

    Since then, the R core team and the community has invested tons of time to improve R and also to make it faster. There are things like RCPP and parallel computing to speed up loops.

    But what still relatively few R users know: loops are not that slow anymore. We want to demonstrate this using two examples.

    Example 1: sqrt()

    We use three ways to calculate the square root of a vector of random numbers:

    1. Vectorized calculation. This will be the way to go because it is internally optimized in C.
    2. A loop. This must be super slow for large vectors.
    3. vapply() (as safe alternative to sapply).

    The three approaches are then compared via bench::mark() regarding their speed for different numbers n of vector lengths. The results are then compared first regarding absolute median times, and secondly (using an independent run), on a relative scale (1 is the vectorized approach). 

    library(tidyverse)
library(bench)

# Calculate square root for each element in loop
sqrt_loop <- function(x) {
  out <- numeric(length(x))
  for (i in seq_along(x)) {
    out[i] <- sqrt(x[i])
  }
  out
}

# Example
sqrt_loop(1:4) # 1.000000 1.414214 1.732051 2.000000

# Compare its performance with two alternatives
sqrt_benchmark <- function(n) {
  x <- rexp(n)
  mark(
    vectorized = sqrt(x),
    loop = sqrt_loop(x),
    vapply = vapply(x, sqrt, FUN.VALUE = 0.0),
    # relative = TRUE
  )
}

# Combine results of multiple benchmarks and plot results
multiple_benchmarks <- function(one_bench, N) {
  res <- vector("list", length(N))
  for (i in seq_along(N)) {
    res[[i]] <- one_bench(N[i]) %>% 
      mutate(n = N[i], expression = names(expression))
  }
  
  ggplot(bind_rows(res), aes(n, median, color = expression)) +
    geom_point(size = 3) +
    geom_line(size = 1) +
    scale_x_log10() +
    ggtitle(deparse1(substitute(one_bench))) +
    theme(legend.position = c(0.8, 0.15))
}

# Apply simulation
multiple_benchmarks(sqrt_benchmark, N = 10^seq(3, 6, 0.25))

    Absolute timings

    Absolute median times on the “sqrt()” task

    Relative timings (using a second run)

    Relative median times of a separate run on the “sqrt()” task

    We see:

    • Run times increase quite linearly with vector size.
    • Vectorization is more than ten times faster than the naive loop.
    • Most strikingly, vapply() is much slower than the naive loop. Would you have thought this?

    Example 2: paste()

    For the second example, we use a less simple function, namely

    paste(“Number”, prettyNum(x, digits = 5))

    What will our three approaches (vectorized, naive loop, vapply) show on this task?

    pretty_paste <- function(x) {
  paste("Number", prettyNum(x, digits = 5))
}

# Example
pretty_paste(pi) # "Number 3.1416"

# Again, call pretty_paste() for each element in a loop
paste_loop <- function(x) {
  out <- character(length(x))
  for (i in seq_along(x)) {
    out[i] <- pretty_paste(x[i])
  }
  out
}

# Compare its performance with two alternatives
paste_benchmark <- function(n) {
  x <- rexp(n)
  mark(
    vectorized = pretty_paste(x),
    loop = paste_loop(x),
    vapply = vapply(x, pretty_paste, FUN.VALUE = ""),
    # relative = TRUE
  )
}

multiple_benchmarks(paste_benchmark, N = 10^seq(3, 5, 0.25))

    Absolute timings

    Absolute median times on the “paste()” task

    Relative timings (using a second run)

    Relative median times of a separate run on the “paste()” task
    • In contrast to the first example, vapply() is now as fast as the naive loop.
    • The time advantage of the vectorized approach is much less impressive. The loop takes in median only 50% longer.

    Conclusion

    1. Vectorization is fast and easy to read. If available, use this. No surprise.
    2. If you use vapply/sapply/lapply, do it for the style, not for the speed. In some cases, the loop will be faster. And, depending on the situation and the audience, a loop might actually be even easier to read.

    The code can be found on github.

    The runs have been made on a Windows 11 system with a four core Intel(R) Core(TM) i7-8650U CPU @ 1.90GHz processor.