• Typical programs operate lines sequentially:

# Import packages
import numpy as np

# Define an array of numbers
foo = np.array([0, 1, 2, 3, 4, 5])

# Define a function that squares numbers
def bar(x):
    return x * x

# Loop over each element and perform an action on it
for element in foo:

        # Print the result of bar
        print(bar(element))

0
1
4
9
16
25

• Say we have a number of input files, like .jpg images, that we want to perform the same actions on, like rotate by 180 degrees and convert to a different format
• We can define a function that takes a file as input and performs these actions, then map it on a list of files

# Import Python Image Library functions
from PIL import Image

• Let’s read an image file and display it:

im = Image.open('./data/kings_cross.jpg')
# Display image
im 

• Parallel computing can be much faster than serial computing
• These problems are essentially independent and share no information between them
• The joblib module makes it simple to run these steps together with a single command
• This workflow is limited to running on a single computer (or compute node) since there is no mechanism to communicate outside

• Install joblib and NumPy if you haven’t done so yet
  • !pip install joblib numpy
• Compare the time of the serial and parallel versions of the following function:

def square(x):
    return x**2

# Create a large array to process
numbers = np.arange(1000000)

# Sequential version
%timeit [square(x) for x in numbers]

• Then try the parallel version:

from joblib import Parallel, delayed

# Parallel version
%timeit Parallel(n_jobs=4)(delayed(square)(x) for x in numbers)

• What did you find? Is the parallel version faster?
• Appendix 01

• Dask arrays are lazy, that is, they do not compute anything until you ask for it
• Lazy functions are useful because they allow the user to build up a computation graph (execution plan) before executing it
• So they optimise the computation and save memory before running anything
• For example, let’s slice the Dask array a to get the first 10 rows of the 6th column

a[:10, 5] # first 10 rows of the 6th column

• We use the .compute() method to compute the result

a[:10, 5].compute()

array([ 0.41269347, -0.99110636, -1.64144012,  0.51300249,  0.61533372,
       -0.40617299,  0.39380648,  2.33899664,  0.2314737 ,  0.60942019])

• Dask Array supports many common NumPy operations including:
• Basic arithmetic and scalar mathematics (+, *, exp, log, etc)
• Reductions along axes (sum(), mean(), std())
• Tensor operations and matrix multiplication (tensordot)
• Array slicing and basic indexing
• Axis reordering and transposition

a.sum().compute()

14.378291513908898

• Let’s compare a similar operation in NumPy:

size = 100000000
np_arr = np.random.random(size)
%timeit np_result = np.sqrt(np_arr) + np.sin(np_arr)

526 ms ± 32.7 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

• And in Dask:

da_arr = da.random.random(size, chunks='auto') 
%timeit da_result = da.sqrt(da_arr) + da.sin(da_arr)

1.15 ms ± 31.6 μs per loop (mean ± std. dev. of 7 runs, 1,000 loops each)

• Since the Dask dataframe we are using has a time series structure, we can use the resample method to aggregate the data by a time period
• Let’s resample the data by 1 hour and calculate the mean of the x and y columns

df[["x", "y"]].resample("1h").mean().head(3)

• We can also use the rolling() method to calculate a rolling mean of the data

df[["x", "y"]].rolling(window="24h").mean().head()

• Install dask if you haven’t done so yet

  • !pip install dask

• Find the right chunk size!

• Create a Dask array with 10 million random numbers (or less if you have memory constraints)

• Vary the chunk size and time the following operation:

• Calculate mean(sqrt(x^2)) on the Dask array

• See the code below for an example with three different chunk sizes. Which one worked best for you? Why do you think that is?

import numpy as np
import dask.array as da

size = 10_000_000

# Dask with SMALL chunks
da_data_small = da.random.random(size, chunks=100_000)  # 100 chunks
%timeit da.sqrt(da_data_small**2).mean().compute()

# Dask with MEDIUM chunks
da_data_medium = da.random.random(size, chunks=2_000_000)  # 5 chunks
%timeit da.sqrt(da_data_medium**2).mean().compute()

# Dask with LARGE chunks
da_data_large = da.random.random(size, chunks=5_000_000)  # Only 2 chunks
%timeit da.sqrt(da_data_large**2).mean().compute()

• Appendix 02

• dask-sql is a library that allows you to run SQL queries on Dask DataFrames
• It is built on top of Apache Calcite, a SQL interpreter
• It is still under development, but it already has many features
• You can install it with pip install dask-sql
• A dask_sql.Context is the Python equivalent to a SQL database,
• It serves as an interface to register all tables and functions used in SQL queries, as well as execute the queries themselves
• A single Context is created and used for the duration of a Python script or notebook

from dask_sql import Context
c = Context()

• Once a Context has been created, there are many ways to register tables in it
• The simplest way to do this is through the create_table method
• Supported input types include Pandas DataFrames, Dask DataFrames, remote datasets (like on Amazon S3), and more
  • More information here: https://dask-sql.readthedocs.io/en/latest/data_input.html

# Register and persist a dask table
from dask.datasets import timeseries
ddf = timeseries()
c.create_table("dask", ddf, persist=True)

• Now we can run SQL queries on the timeseries table

c.sql("SELECT AVG(x) FROM dask").compute()

• You can of course combine all three libraries!

# Perform a multi-column sort
res = c.sql("""
    SELECT * FROM dask ORDER BY name ASC, id DESC, x ASC
""")

# Now do some follow groupby aggregations
res.groupby("name").agg({"x": "sum", "y": "mean"}).compute()

• .csv is very common in data science (and for good reasons)
• Pandas reads and writes .csv files very well, but it is not the best option for large files
• It may need several gigabytes of memory to read a large file, as it reads the entire file into memory
• Dask provides a much more efficient way to read and write .csv files
• Let’s split our dataset

df = dask.datasets.timeseries()
df

import os
import datetime

if not os.path.exists('data'):
    os.mkdir('data')

def name(i):
    """ Provide date for filename given index

    Examples
    --------
    >>> name(0)
    '2000-01-01'
    >>> name(10)
    '2000-01-11'
    """
    return str(datetime.date(2000, 1, 1) + i * datetime.timedelta(days=1))

df.to_csv('data/*.csv', name_function=name);

• We now have many CSV files in our data directory, one for each day in the month of January 2000
• Each CSV file holds time series data for that day
• We can read all of them as one logical dataframe using the dd.read_csv function

df = dd.read_csv('data/2000-*-*.csv')
df

• Let’s do a simple computation

%timeit df.groupby('name').x.mean().compute()

534 ms ± 127 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

• Although .csv files are nice, new formats like Parquet are gaining popularity
• Data are stored by column rather than by row, allowing for efficient querying of specific columns without reading entire rows
• It can be used with various programming languages and data processing frameworks
• Achieves up to 75% reduction in file size compared to uncompressed formats

df.to_parquet('data/2000-01.parquet', engine='pyarrow')

• Now we can read the parquet file

df = dd.read_parquet('data/2000-01.parquet', columns=['name', 'x'], engine='pyarrow')
df

• Let’s do the same computation as before

%timeit df.groupby('name').x.mean().compute()

186 ms ± 11.8 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)

• We can even visualise the computation graph

@dask.delayed
def generate_data(size):
    return np.random.rand(size)

@dask.delayed
def transform_data(data):
    return np.sqrt(data) + np.sin(data)

@dask.delayed
def aggregate_data(data):
    return {
        'mean': np.mean(data),
        'std': np.std(data),
        'max': np.max(data)
    }

# Compare execution
sizes = [1000000, 2000000, 3000000]

# Dask execution
dask_results = []
for size in sizes:
    data = generate_data(size)
    transformed = transform_data(data)
    stats = aggregate_data(transformed)
    dask_results.append(stats)

%timeit dask.compute(*dask_results)

34.2 ms ± 630 μs per loop (mean ± std. dev. of 7 runs, 10 loops each)

dask.visualize(*dask_results)

• Parallel computing can be a powerful tool, but it is not always the best solution
• So here are some tips from the creators of Dask themselves:
• Start small: NumPy, pandas, scikit-learn may have faster functions for what you’re trying to do. Sometimes just by changing the data format to Parquet you can get a huge speedup
• Sampling: If you have a large dataset, try working with a sample of the data first. Do you really need all those TBs of data?
• Load data with Dask: If you have a large dataset, load it with Dask from the start. It’s much easier to scale up than to scale down
• Use .compute() and .persist() sparingly: These functions can be expensive, so use them only when you need to
• Experiment with chunk sizes: The right chunk size can make a big difference in performance (or use chunks='auto' if you’re unsure)
• Break up computations into many pieces: This will allow Dask to parallelise the computation faster
• Use .parquet files for large datasets: They are much more efficient than .csv files
• More tips can be found here

• Here is the solution to the exercise:

def square(x):
    return x**2

# Create a large array to process
numbers = np.arange(1000000)

# Sequential version
%timeit [square(x) for x in numbers]

76.8 ms ± 3.29 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)

from joblib import Parallel, delayed

# Parallel version
%timeit Parallel(n_jobs=4)(delayed(square)(x) for x in numbers)

2.49 s ± 98.5 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

Back to exercise

• Create an array with 10 million random numbers and calculate: mean(sqrt(x^2))
• Try different chunk sizes and see which one works best

size = 10000000

# Dask with SMALL chunks 
da_data_small = da.random.random(size, chunks=100000)
%timeit da.sqrt(da_data_small**2).mean().compute()

35 ms ± 4.24 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)

# Dask with MEDIUM chunks 
da_data_medium = da.random.random(size, chunks=2000000)
%timeit da.sqrt(da_data_medium**2).mean().compute()

20.8 ms ± 2.17 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)

# Dask with LARGE chunks 
da_data_large = da.random.random(size, chunks=5000000)
%timeit da.sqrt(da_data_large**2).mean().compute()

30.4 ms ± 1.3 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)

# Dask with AUTO chunks 
da_data_auto = da.random.random(size, chunks='auto')
%timeit da.sqrt(da_data_auto**2).mean().compute()

51.8 ms ± 2.57 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)

Back to exercise