Analysis of E. coli growth data

(c) 2024 Justin Bois. This work is licensed under a Creative Commons Attribution License CC-BY 4.0. All code contained herein is licensed under an MIT license.

This document was prepared at Caltech with support financial support from the Donna and Benjamin M. Rosen Bioengineering Center.

This tutorial was generated from an Jupyter notebook. You can download the notebook here.

---

import io

import polars as pl
import polars.selectors as cs
import numpy as np
import scipy.optimize

import colorcet
import iqplot

import bokeh.io
import bokeh.plotting
bokeh.io.output_notebook()

Loading BokehJS ...

In this tutorial, we walk through analysis of growth data as acquired in our bulk assay using the plate reader. The analysis presents three main challenges that we have not covered in previous tutorials.

1. We need to wrangle the data exported from the plate reader in a text file to enable easy access and analysis.
2. We need to make plots, some on semilogarithmic scale, to visualize the results.
3. We need to consider portions of the data set in performing regressions, so we need to select and use only parts of the data set.

Data wrangling

Our first step is taking the plate reader data and converting it to an easy format to interpret. As a sample data set, we will use data obtained by the TAs in one of their experiments while developing this module. The data set is in the file 20240422_ecoli_growth_trial_1.txt, which is available here. The first 57 lines of this file are shown below.

Software Version	3.10.06



Experiment File Path:	C:\Users\User\Desktop\20240506_ecoli_full_plate_test.xpt
Protocol File Path:	C:\Users\Public\Documents\Protocols\od600.prt



Plate Number	Plate 1
Date	5/6/2024
Time	8:56:32 AM
Reader Type:	Synergy H1
Reader Serial Number:	2107069
Reading Type	Reader

Procedure Details

Plate Type	96 WELL PLATE (Use plate lid)
Eject plate on completion	
Set Temperature	Setpoint 37�C
Read	Absorbance Endpoint
    Full Plate
    Wavelengths:  600
    Read Speed: Normal,  Delay: 100 msec,  Measurements/Data Point: 8
Start Kinetic	Runtime 24:00:00 (HH:MM:SS), Interval 0:06:00, 241 Reads
    Shake	Linear: Continuous
    Frequency: 1096 cpm (1 mm)
    Read	Absorbance Endpoint
    Full Plate
    Wavelengths:  600
    Read Speed: Normal,  Delay: 100 msec,  Measurements/Data Point: 8
End Kinetic	



Actual Temperature:	37.0

Results
    1	2	3	4	5	6	7	8	9	10	11	12
A	0.091	0.091	0.094	0.098	0.095	0.096	0.086	0.087	0.088	0.099	0.095	0.094	Read 1:600
B	0.103	0.106	0.118	0.087	0.092	0.101	0.095	0.101	0.090	0.097	0.094	0.090	Read 1:600
C	0.093	0.089	0.091	0.088	0.092	0.089	0.087	0.089	0.086	0.089	0.092	0.092	Read 1:600
D	0.096	0.098	0.094	0.089	0.101	0.103	0.093	0.101	0.093	0.096	0.095	0.092	Read 1:600
E	0.089	0.088	0.089	0.085	0.087	0.085	0.095	0.093	0.088	0.094	0.086	0.089	Read 1:600
F	0.098	0.092	0.094	0.092	0.094	0.091	0.100	0.103	0.094	0.097	0.087	0.089	Read 1:600
G	0.085	0.086	0.089	0.092	0.088	0.100	0.089	0.098	0.092	0.099	0.086	0.089	Read 1:600
H	0.087	0.094	0.096	0.088	0.089	0.090	0.087	0.106	0.083	0.099	0.080	0.097	Read 1:600

Read 2:600

Time	T� Read 2:600	A1	A2	A3	A4	A5	A6	A7	A8	A9	A10	A11	A12	B1	B2	B3	B4	B5	B6	B7	B8	B9	B10	B11	B12	C1	C2	C3	C4	C5	C6	C7	C8	C9	C10	C11	C12	D1	D2	D3	D4	D5	D6	D7	D8	D9	D10	D11	D12	E1	E2	E3	E4	E5	E6	E7	E8	E9	E10	E11	E12	F1	F2	F3	F4	F5	F6	F7	F8	F9	F10	F11	F12	G1	G2	G3	G4	G5	G6	G7	G8	G9	G10	G11	G12	H1	H2	H3	H4	H5	H6	H7	H8	H9	H10	H11	H12
0:05:10	37.0	0.090	0.091	0.094	0.105	0.102	0.104	0.088	0.088	0.088	0.097	0.097	0.093	0.112	0.108	0.109	0.092	0.110	0.107	0.104	0.110	0.097	0.095	0.096	0.090	0.093	0.092	0.093	0.093	0.091	0.092	0.090	0.095	0.089	0.095	0.096	0.092	0.102	0.098	0.097	0.097	0.103	0.103	0.100	0.105	0.097	0.098	0.102	0.094	0.094	0.090	0.092	0.092	0.091	0.095	0.094	0.098	0.093	0.094	0.093	0.097	0.098	0.098	0.097	0.097	0.095	0.096	0.093	0.105	0.094	0.097	0.094	0.096	0.094	0.105	0.098	0.096	0.094	0.097	0.095	0.105	0.094	0.102	0.096	0.091	0.096	0.094	0.093	0.092	0.093	0.094	0.089	0.109	0.083	0.098	0.086	0.100
0:11:10	37.0	0.090	0.090	0.094	0.105	0.102	0.104	0.086	0.088	0.088	0.096	0.096	0.093	0.110	0.107	0.109	0.092	0.110	0.106	0.104	0.109	0.097	0.095	0.095	0.090	0.093	0.091	0.092	0.092	0.091	0.091	0.090	0.095	0.089	0.095	0.096	0.093	0.101	0.097	0.096	0.096	0.103	0.102	0.100	0.102	0.096	0.097	0.101	0.094	0.094	0.090	0.091	0.092	0.091	0.095	0.093	0.097	0.092	0.093	0.093	0.095	0.096	0.096	0.096	0.096	0.095	0.095	0.093	0.104	0.094	0.097	0.093	0.095	0.094	0.103	0.098	0.096	0.094	0.097	0.095	0.104	0.094	0.101	0.095	0.091	0.100	0.093	0.093	0.092	0.093	0.093	0.090	0.109	0.084	0.097	0.085	0.099
0:17:10	37.0	0.090	0.090	0.094	0.105	0.103	0.104	0.086	0.088	0.087	0.096	0.096	0.093	0.109	0.107	0.109	0.092	0.110	0.106	0.104	0.109	0.097	0.095	0.095	0.090	0.093	0.091	0.092	0.092	0.091	0.091	0.090	0.095	0.089	0.098	0.096	0.093	0.101	0.097	0.096	0.096	0.102	0.102	0.100	0.102	0.095	0.097	0.101	0.093	0.094	0.089	0.091	0.092	0.091	0.094	0.093	0.097	0.092	0.093	0.093	0.095	0.096	0.096	0.096	0.096	0.095	0.095	0.092	0.104	0.093	0.097	0.093	0.095	0.093	0.104	0.098	0.096	0.094	0.097	0.095	0.104	0.094	0.101	0.090	0.091	0.102	0.093	0.092	0.092	0.093	0.093	0.090	0.110	0.087	0.098	0.086	0.100

Evidently, the plate reader returns a text file with some useful, but unformatted metadata about the measurement, followed by the data. Following the metadata is the first set of absorbance measurements in a table. Following that is the text Read 2:600, which contains a tab-delimited table of absorbance measurements.

Pulling the first time point out is not important for our measurements (the plate reader has not come to temperature anyway), so we will neglect it. Therefore, to pull the data set out of the file, we need to load in the last table. The information we need to do this are the file name, the row number (with zero-based indexing) of the header for the table (containing the well labels), and the last row of the file containing data. We can get all of this information by opening the file in a text editor. For this particular data set, these parameters are defined below.

# Useful info from file from plate reader
fname = "20240514_ecoli_growth_1.txt"
header_row = 53
last_row = 294

We can now write a function to pull the data out and store them in a data frame. Importantly, we will add a column labeled 'time (s)' which is the amount of time, in seconds, since the first measurement we consider, and a column labeled 'time (hr)', which is the same time, but in units of hours.

def timestamp_to_seconds(t_stamp):
    """Convert H:M:S string to seconds. This is necessary
    to hand-build because the datetime utilities do not handle
    times over 24 hours."""
    h, m, s = t_stamp.split(':')
    return 3600 * int(h) + 60 * int(m) + float(s)


def read_plate_reader(fname, header_row, last_row):
    """Read in time series data from Bi 1x plate reader.

    Parameters
    ----------
    fname : str
        Path to file outputted from plate reader.
    header_row : int
        Number of row containing header column (zero-indexed)
    last_row : int
        Last row of the file containing the time series data (zero-indexed)

    Returns
    -------
    output : DataFrame
        Polars DataFrame containing time series data.
    """
    # Get relevant contents of file as string
    with open(fname, "r", encoding="ISO-8859-1") as f:
        lines = f.readlines()
    file_str = "".join(lines[header_row : last_row + 1])

    # How many total lines to read
    n_rows = last_row - header_row

    # Read in data frame
    df = pl.read_csv(
        io.StringIO(file_str),
        separator="\t",
        encoding="ISO-8859-1",
    )

    # Rename temperature column
    T_col = [col for col in df.columns if "T°" in col][0]
    df = df.rename({T_col: "temperature (deg C)"})

    # Parse the time column
    df = df.with_columns(
        pl.col("Time")
        .map_elements(timestamp_to_seconds, return_dtype=float)
        .alias("time (s)")
    ).with_columns(
        pl.col('time (s)') - pl.col('time (s)').min()
    )

    # Time in units of hours
    df = df.with_columns((pl.col('time (s)') / 3600).alias('time (hr)'))

    return df

Let's use this function to load in the data and take a look.

df = read_plate_reader(fname, header_row, last_row)

df

shape: (241, 100)

Time	temperature (deg C)	A1	A2	A3	A4	A5	A6	A7	A8	A9	A10	A11	A12	B1	B2	B3	B4	B5	B6	B7	B8	B9	B10	B11	B12	C1	C2	C3	C4	C5	C6	C7	C8	C9	C10	C11	…	F2	F3	F4	F5	F6	F7	F8	F9	F10	F11	F12	G1	G2	G3	G4	G5	G6	G7	G8	G9	G10	G11	G12	H1	H2	H3	H4	H5	H6	H7	H8	H9	H10	H11	H12	time (s)	time (hr)
str	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	…	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64	f64
"0:05:10"	37.0	0.08	0.083	0.086	0.118	0.122	0.121	0.097	0.089	0.086	0.113	0.117	0.117	0.111	0.116	0.085	0.115	0.125	0.133	0.101	0.113	0.111	0.102	0.112	0.111	0.112	0.13	0.108	0.111	0.107	0.127	0.114	0.115	0.124	0.109	0.125	…	0.107	0.106	0.115	0.111	0.116	0.109	0.113	0.127	0.115	0.11	0.111	0.11	0.107	0.121	0.12	0.11	0.11	0.111	0.118	0.115	0.111	0.107	0.095	0.202	0.214	0.134	0.131	0.108	0.11	0.159	0.239	0.218	0.189	0.115	0.13	0.0	0.0
"0:11:10"	37.1	0.081	0.083	0.086	0.119	0.122	0.122	0.097	0.089	0.086	0.112	0.116	0.117	0.118	0.12	0.086	0.113	0.122	0.143	0.101	0.113	0.117	0.103	0.112	0.11	0.113	0.122	0.11	0.112	0.108	0.124	0.116	0.116	0.126	0.107	0.124	…	0.106	0.105	0.114	0.111	0.118	0.109	0.112	0.121	0.115	0.11	0.111	0.108	0.107	0.12	0.124	0.108	0.106	0.11	0.116	0.123	0.11	0.108	0.098	0.202	0.211	0.133	0.13	0.109	0.111	0.161	0.215	0.231	0.189	0.114	0.119	360.0	0.1
"0:17:10"	37.0	0.083	0.083	0.086	0.121	0.123	0.123	0.097	0.089	0.086	0.113	0.116	0.115	0.101	0.118	0.086	0.114	0.122	0.121	0.102	0.114	0.112	0.104	0.113	0.111	0.113	0.12	0.114	0.112	0.11	0.125	0.12	0.111	0.124	0.112	0.126	…	0.105	0.105	0.114	0.111	0.118	0.109	0.113	0.12	0.116	0.111	0.113	0.11	0.107	0.119	0.12	0.11	0.106	0.111	0.12	0.127	0.104	0.115	0.101	0.202	0.215	0.133	0.13	0.109	0.112	0.164	0.215	0.224	0.187	0.114	0.113	720.0	0.2
"0:23:10"	37.0	0.083	0.083	0.086	0.123	0.126	0.126	0.097	0.089	0.086	0.114	0.118	0.117	0.102	0.117	0.085	0.114	0.123	0.136	0.102	0.114	0.112	0.102	0.115	0.112	0.116	0.117	0.112	0.113	0.111	0.127	0.116	0.113	0.126	0.113	0.126	…	0.106	0.105	0.116	0.112	0.119	0.11	0.113	0.126	0.117	0.113	0.117	0.11	0.107	0.119	0.123	0.108	0.106	0.112	0.131	0.129	0.103	0.111	0.104	0.202	0.209	0.134	0.133	0.111	0.113	0.158	0.214	0.224	0.17	0.114	0.115	1080.0	0.3
"0:29:10"	37.1	0.081	0.082	0.086	0.126	0.129	0.129	0.097	0.089	0.086	0.115	0.119	0.118	0.102	0.141	0.085	0.116	0.125	0.144	0.103	0.116	0.113	0.103	0.115	0.114	0.115	0.118	0.113	0.115	0.113	0.128	0.117	0.123	0.128	0.11	0.128	…	0.105	0.106	0.122	0.113	0.116	0.111	0.115	0.127	0.118	0.113	0.111	0.11	0.107	0.119	0.124	0.11	0.106	0.113	0.133	0.132	0.103	0.112	0.105	0.203	0.211	0.135	0.13	0.112	0.114	0.16	0.215	0.22	0.167	0.114	0.115	1440.0	0.4
…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…	…
"23:41:10"	37.1	1.691	0.083	0.085	1.771	1.774	1.777	0.097	0.088	0.086	0.9	1.121	1.124	0.922	0.913	0.085	1.206	1.236	1.026	0.991	1.116	1.048	1.086	1.132	1.149	1.17	1.09	1.088	1.064	1.059	1.068	1.066	1.087	1.081	1.095	1.116	…	0.919	0.903	1.025	1.035	1.026	0.995	1.024	1.04	0.99	0.984	1.04	0.679	0.618	0.63	1.032	0.942	0.82	1.082	1.155	0.79	0.887	1.09	1.037	0.293	0.299	0.147	0.488	1.056	1.077	0.763	0.963	0.297	0.254	1.065	1.106	84960.0	23.6
"23:47:10"	37.1	1.69	0.082	0.085	1.77	1.773	1.776	0.097	0.088	0.086	0.9	1.119	1.122	0.92	0.913	0.085	1.204	1.228	1.026	0.99	1.114	1.047	1.086	1.131	1.147	1.169	1.089	1.087	1.062	1.057	1.067	1.065	1.087	1.08	1.093	1.115	…	0.918	0.902	1.024	1.035	1.025	0.995	1.023	1.039	0.99	0.983	1.04	0.679	0.618	0.631	1.032	0.941	0.825	1.082	1.155	0.79	0.887	1.089	1.035	0.293	0.299	0.145	0.488	1.055	1.075	0.762	0.964	0.282	0.254	1.064	1.105	85320.0	23.7
"23:53:10"	37.0	1.688	0.083	0.085	1.77	1.773	1.775	0.097	0.088	0.086	0.899	1.118	1.121	0.92	0.912	0.085	1.203	1.228	1.036	0.99	1.114	1.046	1.087	1.13	1.147	1.17	1.089	1.086	1.061	1.056	1.066	1.064	1.085	1.079	1.092	1.115	…	0.917	0.902	1.025	1.035	1.025	0.995	1.024	1.04	0.989	0.984	1.04	0.678	0.618	0.63	1.032	0.942	0.819	1.082	1.155	0.79	0.888	1.093	1.036	0.324	0.3	0.144	0.485	1.054	1.076	0.764	0.963	0.3	0.254	1.064	1.109	85680.0	23.8
"23:59:10"	37.0	1.687	0.083	0.085	1.769	1.772	1.774	0.097	0.088	0.086	0.899	1.119	1.122	0.919	0.912	0.085	1.204	1.233	1.025	0.99	1.114	1.045	1.087	1.13	1.147	1.168	1.087	1.085	1.06	1.057	1.065	1.063	1.084	1.079	1.092	1.114	…	0.916	0.901	1.026	1.034	1.025	0.995	1.023	1.039	0.988	0.983	1.039	0.678	0.617	0.63	1.033	0.943	0.82	1.082	1.155	0.79	0.889	1.094	1.036	0.294	0.299	0.145	0.487	1.054	1.076	0.762	0.963	0.296	0.254	1.063	1.104	86040.0	23.9
"24:05:10"	37.0	1.687	0.083	0.085	1.768	1.771	1.773	0.097	0.089	0.086	0.898	1.118	1.121	0.918	0.912	0.085	1.203	1.226	1.024	0.99	1.113	1.045	1.087	1.129	1.146	1.168	1.087	1.084	1.06	1.056	1.064	1.062	1.084	1.078	1.091	1.114	…	0.915	0.899	1.025	1.034	1.025	0.994	1.022	1.038	0.988	0.982	1.046	0.676	0.617	0.63	1.032	0.942	0.82	1.083	1.156	0.789	0.888	1.09	1.038	0.293	0.301	0.144	0.485	1.053	1.074	0.761	0.964	0.3	0.254	1.062	1.104	86400.0	24.0

We now conveniently have a data set with our results.

We are, however, missing some metadata. We want to know, e.g., what conditions were present in well C6, for example. We can create a dictionary for our samples.

well_metadata = dict(
    A1="LB blank",
    A2="LB blank",
    A3="LB blank",
    A4="LB",
    A5="LB",
    A6="LB",
    A7="M9 blank",
    A8="M9 blank",
    A9="M9 blank",
    A10="glucose",
    A11="glucose",
    A12="glucose",
    B1="xylose",
    B2="xylose",
    B3="xylose",
    B4="xylose:glucose 4:1",
    B5="xylose:glucose 4:1",
    B6="xylose:glucose 4:1",
    B7="xylose:glucose 2:1",
    B8="xylose:glucose 2:1",
    B9="xylose:glucose 2:1",
    B10="xylose:glucose 1:1",
    B11="xylose:glucose 1:1",
    B12="xylose:glucose 1:1",
    C1="sorbitol",
    C2="sorbitol",
    C3="sorbitol",
    C4="sorbitol:glucose 4:1",
    C5="sorbitol:glucose 4:1",
    C6="sorbitol:glucose 4:1",
    C7="sorbitol:glucose 2:1",
    C8="sorbitol:glucose 2:1",
    C9="sorbitol:glucose 2:1",
    C10="sorbitol:glucose 1:1",
    C11="sorbitol:glucose 1:1",
    C12="sorbitol:glucose 1:1",
    D1="galactose",
    D2="galactose",
    D3="galactose",
    D4="galactose:glucose 4:1",
    D5="galactose:glucose 4:1",
    D6="galactose:glucose 4:1",
    D7="galactose:glucose 2:1",
    D8="galactose:glucose 2:1",
    D9="galactose:glucose 2:1",
    D10="galactose:glucose 1:1",
    D11="galactose:glucose 1:1",
    D12="galactose:glucose 1:1",
    E1="rhamnose",
    E2="rhamnose",
    E3="rhamnose",
    E4="rhamnose:glucose 4:1",
    E5="rhamnose:glucose 4:1",
    E6="rhamnose:glucose 4:1",
    E7="rhamnose:glucose 2:1",
    E8="rhamnose:glucose 2:1",
    E9="rhamnose:glucose 2:1",
    E10="rhamnose:glucose 1:1",
    E11="rhamnose:glucose 1:1",
    E12="rhamnose:glucose 1:1",
    F1="ribose",
    F2="ribose",
    F3="ribose",
    F4="ribose:glucose 4:1",
    F5="ribose:glucose 4:1",
    F6="ribose:glucose 4:1",
    F7="ribose:glucose 2:1",
    F8="ribose:glucose 2:1",
    F9="ribose:glucose 2:1",
    F10="ribose:glucose 1:1",
    F11="ribose:glucose 1:1",
    F12="ribose:glucose 1:1",
    G1="arabinose",
    G2="arabinose",
    G3="arabinose",
    G4="arabinose:glucose 4:1",
    G5="arabinose:glucose 4:1",
    G6="arabinose:glucose 4:1",
    G7="arabinose:glucose 2:1",
    G8="arabinose:glucose 2:1",
    G9="arabinose:glucose 2:1",
    G10="arabinose:glucose 1:1",
    G11="arabinose:glucose 1:1",
    G12="arabinose:glucose 1:1",
    H1="student choice 1",
    H2="student choice 2",
    H3="student choice 3",
    H4="student choice 4",
    H5="student choice 5",
    H6="student choice 6",
    H7="student choice 7",
    H8="student choice 8",
    H9="student choice 9",
    H10="student choice 10",
    H11="student choice 11",
    H12="student choice 12",
)

It is convenient to invert the dictionary so that we can look at which wells correspond to a given condition (e.g., wells B2 and C2 have glucose).

def invert_dict(input_dict):
    """Invert a dictionary ignore None values."""
    # Create an empty dictionary to store the inverted mapping
    output_dict = {}

    # Iterate over items in the input dictionary
    for key, value in input_dict.items():
        if value is not None:
            if value in output_dict:
                output_dict[value].append(key)
            else:
                output_dict[value] = [key]

    return output_dict

Let's now make and check out our inverse dictionary.

conditions = invert_dict(well_metadata)

conditions

{'LB blank': ['A1', 'A2', 'A3'],
 'LB': ['A4', 'A5', 'A6'],
 'M9 blank': ['A7', 'A8', 'A9'],
 'glucose': ['A10', 'A11', 'A12'],
 'xylose': ['B1', 'B2', 'B3'],
 'xylose:glucose 4:1': ['B4', 'B5', 'B6'],
 'xylose:glucose 2:1': ['B7', 'B8', 'B9'],
 'xylose:glucose 1:1': ['B10', 'B11', 'B12'],
 'sorbitol': ['C1', 'C2', 'C3'],
 'sorbitol:glucose 4:1': ['C4', 'C5', 'C6'],
 'sorbitol:glucose 2:1': ['C7', 'C8', 'C9'],
 'sorbitol:glucose 1:1': ['C10', 'C11', 'C12'],
 'galactose': ['D1', 'D2', 'D3'],
 'galactose:glucose 4:1': ['D4', 'D5', 'D6'],
 'galactose:glucose 2:1': ['D7', 'D8', 'D9'],
 'galactose:glucose 1:1': ['D10', 'D11', 'D12'],
 'rhamnose': ['E1', 'E2', 'E3'],
 'rhamnose:glucose 4:1': ['E4', 'E5', 'E6'],
 'rhamnose:glucose 2:1': ['E7', 'E8', 'E9'],
 'rhamnose:glucose 1:1': ['E10', 'E11', 'E12'],
 'ribose': ['F1', 'F2', 'F3'],
 'ribose:glucose 4:1': ['F4', 'F5', 'F6'],
 'ribose:glucose 2:1': ['F7', 'F8', 'F9'],
 'ribose:glucose 1:1': ['F10', 'F11', 'F12'],
 'arabinose': ['G1', 'G2', 'G3'],
 'arabinose:glucose 4:1': ['G4', 'G5', 'G6'],
 'arabinose:glucose 2:1': ['G7', 'G8', 'G9'],
 'arabinose:glucose 1:1': ['G10', 'G11', 'G12'],
 'student choice 1': ['H1'],
 'student choice 2': ['H2'],
 'student choice 3': ['H3'],
 'student choice 4': ['H4'],
 'student choice 5': ['H5'],
 'student choice 6': ['H6'],
 'student choice 7': ['H7'],
 'student choice 8': ['H8'],
 'student choice 9': ['H9'],
 'student choice 10': ['H10'],
 'student choice 11': ['H11'],
 'student choice 12': ['H12']}

These kind of dictionaries need to be hand-built, based on what you do as an experimenter. For convenience, below is code for making a dictionary for your experiment. You should update the wells in row H to reflect what you and your classmates decided to put in there.

well_metadata = dict(
    A1="LB blank",
    A2="LB blank",
    A3="LB blank",
    A4="LB",
    A5="LB",
    A6="LB",
    A7="M9 blank",
    A8="M9 blank",
    A9="M9 blank",
    A10="glucose",
    A11="glucose",
    A12="glucose",
    B1="xylose",
    B2="xylose",
    B3="xylose",
    B4="xylose:glucose 4:1",
    B5="xylose:glucose 4:1",
    B6="xylose:glucose 4:1",
    B7="xylose:glucose 2:1",
    B8="xylose:glucose 2:1",
    B9="xylose:glucose 2:1",
    B10="xylose:glucose 1:1",
    B11="xylose:glucose 1:1",
    B12="xylose:glucose 1:1",
    C1="sorbitol",
    C2="sorbitol",
    C3="sorbitol",
    C4="sorbitol:glucose 4:1",
    C5="sorbitol:glucose 4:1",
    C6="sorbitol:glucose 4:1",
    C7="sorbitol:glucose 2:1",
    C8="sorbitol:glucose 2:1",
    C9="sorbitol:glucose 2:1",
    C10="sorbitol:glucose 1:1",
    C11="sorbitol:glucose 1:1",
    C12="sorbitol:glucose 1:1",
    D1="galactose",
    D2="galactose",
    D3="galactose",
    D4="galactose:glucose 4:1",
    D5="galactose:glucose 4:1",
    D6="galactose:glucose 4:1",
    D7="galactose:glucose 2:1",
    D8="galactose:glucose 2:1",
    D9="galactose:glucose 2:1",
    D10="galactose:glucose 1:1",
    D11="galactose:glucose 1:1",
    D12="galactose:glucose 1:1",
    E1="rhamnose",
    E2="rhamnose",
    E3="rhamnose",
    E4="rhamnose:glucose 4:1",
    E5="rhamnose:glucose 4:1",
    E6="rhamnose:glucose 4:1",
    E7="rhamnose:glucose 2:1",
    E8="rhamnose:glucose 2:1",
    E9="rhamnose:glucose 2:1",
    E10="rhamnose:glucose 1:1",
    E11="rhamnose:glucose 1:1",
    E12="rhamnose:glucose 1:1",
    F1="ribose",
    F2="ribose",
    F3="ribose",
    F4="ribose:glucose 4:1",
    F5="ribose:glucose 4:1",
    F6="ribose:glucose 4:1",
    F7="ribose:glucose 2:1",
    F8="ribose:glucose 2:1",
    F9="ribose:glucose 2:1",
    F10="ribose:glucose 1:1",
    F11="ribose:glucose 1:1",
    F12="ribose:glucose 1:1",
    G1="lactose",
    G2="lactose",
    G3="lactose",
    G4="lactose:glucose 4:1",
    G5="lactose:glucose 4:1",
    G6="lactose:glucose 4:1",
    G7="lactose:glucose 2:1",
    G8="lactose:glucose 2:1",
    G9="lactose:glucose 2:1",
    G10="lactose:glucose 1:1",
    G11="lactose:glucose 1:1",
    G12="lactose:glucose 1:1",
    H1="student choice 1",
    H2="student choice 2",
    H3="student choice 3",
    H4="student choice 4",
    H5="student choice 5",
    H6="student choice 6",
    H7="student choice 7",
    H8="student choice 8",
    H9="student choice 9",
    H10="student choice 10",
    H11="student choice 11",
    H12="student choice 12",
)

Plotting data

We now have a nice, tidy data set to work with. Let's start making some plots.

Temperature check

First, we'll plot the temperature versus time to make sure there we no problems with the temperature setting in the plate reader.

p = bokeh.plotting.figure(
    frame_width=400,
    frame_height=200,
    x_axis_label='time (hr)',
    y_axis_label='temperature (deg. C)',
    x_range=[0, 24],
)

p.line(x=df['time (hr)'], y=df['temperature (deg C)'], line_width=2)

bokeh.io.show(p)

The temperature was held steady, so we're all good there.

Checking blanks and computing the background

Next, let's look at the background traces. We will plot all time series for the blanks, and then overlay the average blank,which we will store as the background.

# Background traces
p = bokeh.plotting.figure(
    frame_width=400,
    frame_height=200,
    x_axis_label='time (hr)',
    y_axis_label='absorbance',
    x_range=[0, 24],
)

# Plot all background traces
for well in conditions['M9 blank']:
    p.line(x=df['time (hr)'], y=df[well])

# Compute background by averaging background values at each time point
df = df.with_columns(pl.mean_horizontal(conditions['M9 blank']).alias('background'))

# Overlay average blank
p.line(x=df['time (hr)'], y=df['background'], color='tomato', line_width=2)

bokeh.io.show(bokeh.layouts.column(p))

All of the blanks are more or less steady, with a typical absorbance of 0.088.

Plotting growth curves in pure sugars

Now, we can begin plotting the background-subtracted absorbance. Let us start with a simple example, pure glocose.

# Instantiate figure
p = bokeh.plotting.figure(
    frame_width=400,
    frame_height=200,
    x_axis_label='time (hr)',
    y_axis_label='absorbance',
    x_range=[0, 24],
)

# Loop through and add a line for each pure glucose trace
for well in conditions['glucose']:
    p.line(x=df['time (hr)'], y=df[well] - df['background'])

bokeh.io.show(p)

This is interesting. We immediately see exponential growth, which ends at about the seven hour mark, when the bacterial abruptly stop growing. There is a pause for a couple hours, another growth phase for about five hours, and then a slow decline. We will postulate in class why there may be an extra growth phase.

We would like to compare with the other pure sugar sources as well, so let's include galactose, sorbitol, and xylose in our plot. We will color galactose red, sorbitol orange, and xylose purple.

# Loop through and add a line for each pure galactose trace
for well in conditions['galactose']:
    p.line(x=df['time (hr)'], y=df[well] - df['background'], color='orange')
    
# Loop through and add a line for each pure sorbitol trace
for well in conditions['sorbitol']:
    p.line(x=df['time (hr)'], y=df[well] - df['background'], color='green')

# Loop through and add a line for each pure xylose trace
for well in conditions['xylose']:
    p.line(x=df['time (hr)'], y=df[well] - df['background'], color='tomato')

bokeh.io.show(p)

All curves appear to show exponential growth. Glucose, xylose, and sorbitol all then have a pause, and then secondary growth and then decay Though note that one of the xylose curves never grew.

This plot is fine and good, but it leaves a bit to be desired. First, we would like to have the absorbance axis on a logarithmic scale so we can better see the exponential growth; it will appear as a line, with the growth rate being the slope. Second, we would like to have a legend, so we do not have to just remember which curve was which. To get a logarithmic vertical axis, we use the y_axis_type='log' keyword argument when we instantiate the plot. There are many ways to get a legend; I prefer hand-building it and then placing it outside the plot. We can also make the legend clickable, so that we can show and hide specific curves. The code below accomplishes all of this.

# Instantiate figure
p = bokeh.plotting.figure(
    frame_width=400,
    frame_height=200,
    x_axis_label="time (hr)",
    y_axis_label="absorbance",
    y_axis_type="log",
    x_range=[0, 24],
    toolbar_location="above",
)

# Items we will place in the legend
legend_items = []

# Loop through conditions with pretty colors
for condition, color in zip(
    ["glucose", "galactose", "sorbitol", "xylose", "ribose", "rhamnose", "arabinose"],
    colorcet.b_glasbey_category10,
):
    # Instantiate plotted lines for a given sugar
    lines = []

    # Loop through each well corresponding to the sugar
    for well in conditions[condition]:
        # A line of the time series of absorbance
        lines.append(
            p.line(x=df["time (hr)"], y=df[well] - df["background"], line_color=color)
        )

    # Add the lines to the legend
    legend_items.append((condition, lines))

# Create the legend from the items
legend = bokeh.models.Legend(items=legend_items, click_policy="hide")

# Add the legend to the plot
p.add_layout(legend, "right")

bokeh.io.show(p)

This plot is clear. There is about a 2.5-hour lag phase for glucose before exponential growth kicks in. The lag is a bit longer for sorbitol and xylose. Sorbitol and xylose seem to have bacteria with similar growth rates, which are smaller than that in glucose. We will quantify these more precisely when we do regressions in a moment.

Plotting multiple sugar sources

Let us now look at plots of growth curves in the presence of two sugars. We will look at glucose and sorbitol to start with.

It is good practice to put all of the curves on the same plot so that they can be easily compared. Since we want to do this for multiple sugar sources, we can write a function to generate the plots.

def plot_sugar_source(df, sugar_source, y_axis_type='log'):
    """Make a plot of all glucose + sugar_source conditions"""
    # Instantiate figure
    p = bokeh.plotting.figure(
        frame_width=450,
        frame_height=250,
        x_axis_label="time (hr)",
        y_axis_label="absorbance",
        y_axis_type=y_axis_type,
        x_range=[0, 24],
        toolbar_location="above",
    )

    # Conditions to consider
    conds = ["glucose", sugar_source] + [
        f"{sugar_source}:glucose {ratio}" for ratio in ["1:1", "2:1", "4:1"]
    ]

    # Color glucose in blue, pure other sugar in orange, grays for ratios
    colors = ["#1f77b4", "orange"] + list(bokeh.palettes.Greys5[1:-1])

    # Populate glyphs, legend items
    legend_items = []
    for cond, color in zip(conds, colors):
        lines = []
        for well in conditions[cond]:
            lines.append(
                p.line(
                    x=df["time (hr)"],
                    y=df[well] - df["background"],
                    line_color=color,
                )
            )
        legend_items.append((cond, lines))

    # Build and add legend
    legend = bokeh.models.Legend(items=legend_items, click_policy="hide")
    p.add_layout(legend, "right")

    return p

Let's take this for a spin! We'll use sorbitol.

bokeh.io.show(plot_sugar_source(df, 'sorbitol'))

This plot clearly shows diauxic growth! All curves follow glucose at short times, and then when the glucose is eaten, they switch to sorbitol with a different slope.

Obtaining the growth rates

We learned how to do a nonlinear regression in a previous tutorial. Here, we need to extract the region of the growth curve of interest (exponential growth phase) and then directly what we applied in the previous tutorial.

Extracting data of interest

To extract data of interest from a data frame, we can use Boolean indexing. Let's look at the pure glucose curves, shown again here.

# Instantiate figure
p = bokeh.plotting.figure(
    frame_width=400,
    frame_height=200,
    x_axis_label='time (hr)',
    y_axis_label='absorbance',
    y_axis_type='log',
    x_range=[0, 24],
)

# Loop through and add a line for each pure glucose trace
for well in conditions['glucose']:
    p.line(x=df['time (hr)'], y=df[well] - df['background'])

bokeh.io.show(p)

In looking at the graph, exponential growth phase is between when the background-subtracted OD is between about 0.1, and 0.4. To extract those parts of the curve for a given well, we use the following syntax. Be sure to note the use of parentheses and the & operator.

OD_range = [0.1, 0.4]
col = conditions['glucose'][0]

sub_df = (
    df
    .filter(
        (pl.col(col) - pl.col('background') >= OD_range[0]) 
        & (pl.col(col) - pl.col('background') <= OD_range[1])
    )
    .select('time (hr)', 'background', col)
)

Now, we have a sub-data frame that only has time points between three and five hours.

Performing a regression to get the growth rate

Let us now consider one of the glucose curves and perform a curve fit. As in the first tutorial, we need to define a theoretical function for exponential growth. Since we are considering background-substracted growth, the equation is

\begin{align} A = A_0\,\mathrm{e}^{r t}, \end{align}

where $A$ is absorbance and $r$ is the growth rate.

# Define exponential growth function
def exp_growth(t, A0, r):
    return A0 * np.exp(r * t)

Next, get guess parameters.

# Parameter guesses
A0_guess = 0.01
r_guess = 1.0 # 1/hr

Next, we pull out the time points and absorbance measurements, making sure to do background subtraction.

# Pull out time points and absorbances
t = sub_df['time (hr)'].to_numpy()
A = (sub_df[conditions['glucose'][0]] - sub_df['background']).to_numpy()

Finally, we use scipy.optimize.curve_fit() to get the best-fit parameters.

# Perform optimization
popt, _ = scipy.optimize.curve_fit(exp_growth, t, A, p0=(A0_guess, r_guess))

# Extract parameters
A0, r = popt

# What was our growth rate?
print(f"growth rate = {r} inverse hours")

growth rate = 0.5933536217655401 inverse hours

A growth rate of 0.6 inverse hours corresponds to a division time of $\ln 2 / 0.6 \approx 1$ hour.

Let's overlay a plot.

# Generator theoretical curve, going a bit beyond
# both sides of curve fit region
t_theor = np.linspace(t.min(), t.max(), 200)
A_theor = exp_growth(t_theor, *popt)

# Add to plot
p.line(t_theor, A_theor, line_color='tomato')

bokeh.io.show(p)

Looks nice!

Growth rates for single sugars

We should automate the regression by writing a function to perform it.

def fit_growth(df, well, OD_range):
    """Obtain estimates for A0 and r for a given sub-trace"""
    # Parameter guesses
    A0_guess = 0.01
    r_guess = 1.0  # 1/hr

    # Pull out data
    sub_df = (
        df
        .filter(
            (pl.col(well) - pl.col('background') >= OD_range[0]) 
            & (pl.col(well) - pl.col('background') <= OD_range[1])
        )
        .select('time (hr)', 'background', well)
    )

    t = sub_df["time (hr)"].to_numpy()
    A = (sub_df[well] - sub_df["background"]).to_numpy()

    # Perform optimization
    popt, _ = scipy.optimize.curve_fit(exp_growth, t, A, p0=(A0_guess, r_guess))

    # Return best fit parameters and also the time range for the fit
    return *popt, [t.min(), t.max()]

We can now use this function to get growth rates for all of the pure sugars. To do so, we will need to specify the OD ranges for the plots above. In your analysis, you will need to obtain these bounds, but for now, we will assume they are all the same, [0.1, 0.4] for a single sugar and [0.1, 0.25], [0.4, 0.6] for diauxic growth.

OD_range = {
    cond: [
        [(0.1, 0.25), (0.4, 0.6)] 
        if ':' in cond 
        else [(0.1, 0.4)]
        for _ in range(len(conditions[cond]))
        ]
    for cond in conditions
}

# Take a look
OD_range

{'LB blank': [[(0.1, 0.4)], [(0.1, 0.4)], [(0.1, 0.4)]],
 'LB': [[(0.1, 0.4)], [(0.1, 0.4)], [(0.1, 0.4)]],
 'M9 blank': [[(0.1, 0.4)], [(0.1, 0.4)], [(0.1, 0.4)]],
 'glucose': [[(0.1, 0.4)], [(0.1, 0.4)], [(0.1, 0.4)]],
 'xylose': [[(0.1, 0.4)], [(0.1, 0.4)], [(0.1, 0.4)]],
 'xylose:glucose 4:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'xylose:glucose 2:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'xylose:glucose 1:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'sorbitol': [[(0.1, 0.4)], [(0.1, 0.4)], [(0.1, 0.4)]],
 'sorbitol:glucose 4:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'sorbitol:glucose 2:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'sorbitol:glucose 1:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'galactose': [[(0.1, 0.4)], [(0.1, 0.4)], [(0.1, 0.4)]],
 'galactose:glucose 4:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'galactose:glucose 2:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'galactose:glucose 1:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'rhamnose': [[(0.1, 0.4)], [(0.1, 0.4)], [(0.1, 0.4)]],
 'rhamnose:glucose 4:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'rhamnose:glucose 2:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'rhamnose:glucose 1:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'ribose': [[(0.1, 0.4)], [(0.1, 0.4)], [(0.1, 0.4)]],
 'ribose:glucose 4:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'ribose:glucose 2:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'ribose:glucose 1:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'arabinose': [[(0.1, 0.4)], [(0.1, 0.4)], [(0.1, 0.4)]],
 'arabinose:glucose 4:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'arabinose:glucose 2:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'arabinose:glucose 1:1': [[(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)],
  [(0.1, 0.25), (0.4, 0.6)]],
 'student choice 1': [[(0.1, 0.4)]],
 'student choice 2': [[(0.1, 0.4)]],
 'student choice 3': [[(0.1, 0.4)]],
 'student choice 4': [[(0.1, 0.4)]],
 'student choice 5': [[(0.1, 0.4)]],
 'student choice 6': [[(0.1, 0.4)]],
 'student choice 7': [[(0.1, 0.4)]],
 'student choice 8': [[(0.1, 0.4)]],
 'student choice 9': [[(0.1, 0.4)]],
 'student choice 10': [[(0.1, 0.4)]],
 'student choice 11': [[(0.1, 0.4)]],
 'student choice 12': [[(0.1, 0.4)]]}

Of course, you'll have to fill out this dictionary with OD ranges that make sense for your analysis. Since we'll play with sorbitol, I'll add good ranges for that here.

OD_range['sorbitol:glucose 4:1'] = [
    [(0.1, 0.2), (0.4, 0.6)],
    [(0.1, 0.2), (0.4, 0.6)],
    [(0.1, 0.2), (0.4, 0.6)]
]

OD_range['sorbitol:glucose 2:1'] = [
    [(0.1, 0.25), (0.45, 0.7)],
    [(0.1, 0.25), (0.45, 0.7)],
    [(0.1, 0.25), (0.45, 0.7)]
]

OD_range['sorbitol:glucose 1:1'] = [
    [(0.1, 0.3), (0.65, 0.8)],
    [(0.1, 0.3), (0.65, 0.8)],
    [(0.1, 0.3), (0.65, 0.8)]
]

OD_range['sorbitol'] = [
    [(0.1, 0.6)],
    [(0.1, 0.6)],
    [(0.1, 0.6)],
]

We will perform the regressions to get the growth rates. We will store each result as a dictionary containing the sugar source, trial number, and value of the growth rate from the regression. Then, we can convert a list of these results dictionaries into a data frame for easy viewing and plotting.

results = []
for sugar in ['glucose', 'galactose', 'sorbitol', 'xylose']:
    for trial, well in enumerate(conditions[sugar]):
        for ODr in OD_range[sugar][trial]:
            try:
                A0, r, _ = fit_growth(df, well, ODr)
                results.append(dict(sugar=sugar, trial=trial, OD_low=ODr[0], OD_high=ODr[1], A0=A0, r=r))
            except:
                print(f'WARNING: {sugar}, {well}, {ODr} has no data in OD range; skipping.')

# Make a data frame from the list of results
df_growth_rates = pl.DataFrame(results)

# Take a look
df_growth_rates

WARNING: xylose, B3, (0.1, 0.4) has no data in OD range; skipping.

shape: (11, 6)

sugar	trial	OD_low	OD_high	A0	r
str	i64	f64	f64	f64	f64
"glucose"	0	0.1	0.4	0.024282	0.593354
"glucose"	1	0.1	0.4	0.025689	0.60446
"glucose"	2	0.1	0.4	0.02751	0.594018
"galactose"	0	0.1	0.4	0.027992	0.229443
"galactose"	1	0.1	0.4	0.026537	0.253469
…	…	…	…	…	…
"sorbitol"	0	0.1	0.6	0.027545	0.392449
"sorbitol"	1	0.1	0.6	0.026191	0.400065
"sorbitol"	2	0.1	0.6	0.024077	0.412325
"xylose"	0	0.1	0.4	0.014733	0.442673
"xylose"	1	0.1	0.4	0.010956	0.449941

If we like, we can summarize our results in a plot.

p = iqplot.strip(
    df_growth_rates,
    q="r",
    cats="sugar",
    q_axis='y',
    y_axis_label="growth rate (1/hr)",
    spread="swarm",
    y_range=[0, 1],
)

bokeh.io.show(p)

This provides a nice look into the growth rates.

As a check, we should make sure the curve fits look ok. We can again overlay them on the plot. So, let's make a function that performs the regressions for a set of conditions and, if desired, makes a plot with the growth curves with regressions overlaid. It's a fair amount of code to do it, but it's straightforward. We will write it as a function, since we will use this again and again.

def multi_reg_plot(
    df,
    conds,
    colors=colorcet.b_glasbey_category10,
    OD_range=None,
    show_regression=True,
    return_plot=True,
):
    """Perform regression for multiple conditions and, if asked,
    generate a plot with multiple conditions, with or without
    showing regression lines."""
    if return_plot:
        # Instantiate figure
        p = bokeh.plotting.figure(
            frame_width=400,
            frame_height=200,
            x_axis_label="time (hr)",
            y_axis_label="absorbance",
            y_axis_type="log",
            x_range=[0, 24],
            toolbar_location="above",
        )
    
        # Items we will place in the legend
        legend_items = []

        # Loop through conditions with pretty colors
        for cond, color in zip(conds, colors):
            # Instantiate plotted lines for a given condition
            lines = []
    
            # Loop through each well corresponding to the sugar
            for well in conditions[cond]:
                # A line of the time series of absorbance
                lines.append(
                    p.line(
                        x=df["time (hr)"], y=df[well] - df["background"], line_color=color
                    )
                )
    
            # Add the lines to the legend
            legend_items.append((cond, lines))

    # Loop through and perform regressions
    fit_lines = []
    results = []
    for cond in conds:
        # Loop through each well corresponding to the sugar
        for trial, well in enumerate(conditions[cond]):
            # Get best fit curve for each OD range
            for diaux_phase, ODr in enumerate(OD_range[cond][trial]):
                try:
                    # Get best fit parameters
                    A0, r, t_range = fit_growth(df, well, ODr)

                    # Time points for theoretical curve
                    t_theor = np.linspace(*t_range, 200)

                    # Theoretical absorbance
                    A_theor = exp_growth(t_theor, A0, r)

                    results.append(
                        dict(
                            condition=cond,
                            trial=trial,
                            diaux_phase=diaux_phase,
                            OD_low=ODr[0],
                            OD_high=ODr[1],
                            A0=A0,
                            r=r,
                        )
                    )

                    # Plot theoretical curve in gray
                    if return_plot and show_regression:
                        fit_lines.append(
                            p.line(x=t_theor, y=A_theor, line_color="gray")
                        )
                except:
                    print(
                        f"WARNING: {sugar}, {well}, {ODr} has no data in OD range; skipping."
                    )

    # Add the fit lines to the legend
    if return_plot and show_regression:
        legend_items.append(("regression", fit_lines))

    # Make a data frame from the list of results
    df_growth_rates = pl.DataFrame(results)

    if return_plot:
        # Create the legend from the items
        legend = bokeh.models.Legend(items=legend_items, click_policy="hide")
    
        # Add the legend to the plot
        p.add_layout(legend, "right")

        return df_growth_rates, p
    else:
        return p

Let's now do our regressions and make our plot to do the check!

df_pure_sugars, p = multi_reg_plot(
    df,
    ["glucose", "galactose", "sorbitol", "xylose"],
    OD_range=OD_range,
)

bokeh.io.show(p)

WARNING: xylose, B3, (0.1, 0.4) has no data in OD range; skipping.

Fitting diauxic curves

Fitting a diauxic curve is a bit more challenging, since we need to look carefully at the respective exponential growth regimes. Let's consider again the sorbitol curve as an example.

df_sorb, p = multi_reg_plot(
    df,
    ["glucose", "sorbitol", 'sorbitol:glucose 4:1', 'sorbitol:glucose 2:1', 'sorbitol:glucose 1:1'],
    colors=['#1f77b4', 'orange'] + list(bokeh.palettes.Purples5[1:-1]),
    OD_range=OD_range,
)


bokeh.io.show(p)

The regressions look good! Now let's make a plot of the growth rates we got from the regression.

# Look at result
df_sorb

shape: (24, 7)

condition	trial	diaux_phase	OD_low	OD_high	A0	r
str	i64	i64	f64	f64	f64	f64
"glucose"	0	0	0.1	0.4	0.024282	0.593354
"glucose"	1	0	0.1	0.4	0.025689	0.60446
"glucose"	2	0	0.1	0.4	0.02751	0.594018
"sorbitol"	0	0	0.1	0.6	0.027545	0.392449
"sorbitol"	1	0	0.1	0.6	0.026191	0.400065
…	…	…	…	…	…	…
"sorbitol:glucose 1:1"	0	1	0.65	0.8	0.101014	0.281996
"sorbitol:glucose 1:1"	1	0	0.1	0.3	0.026459	0.597242
"sorbitol:glucose 1:1"	1	1	0.65	0.8	0.101654	0.298754
"sorbitol:glucose 1:1"	2	0	0.1	0.3	0.026599	0.561807
"sorbitol:glucose 1:1"	2	1	0.65	0.8	0.105969	0.278731

Looks good; now let's make a plot!

# Specify the order of the conditions
order = (
    ("glucose", 0),
    ("sorbitol", 0),
    ("sorbitol:glucose 4:1", 0),
    ("sorbitol:glucose 4:1", 1),
    ("sorbitol:glucose 2:1", 0),
    ("sorbitol:glucose 2:1", 1),
    ("sorbitol:glucose 1:1", 0),
    ("sorbitol:glucose 1:1", 1),
)

# Colors
colors = ['#1f77b4', 'orange'] * 4

# Make a strip plot of the growth rates
p = iqplot.strip(
    df_sorb,
    q="r",
    cats=["condition", "diaux_phase"],
    y_axis_label="growth rate (1/hr)",
    spread="swarm",
    q_axis="y",
    order=order,
    palette=colors,
    frame_width=500,
    y_range=[0, 0.8],
)

bokeh.io.show(p)

We see that, indeed, all cells seem to be growing on glucose at the same rate as in pure glucose, but as glucose runs out, the cells grow roughly at the same rate as in pure sorbitol, though with a decreasing rate as the sorbitol:glucose ratio drops.

Further analysis

It is up to you to look carefully at the growth curves and set the appropriate OD ranges for the regressions and do analysis for the rest of the sugars and diauxic combinations, as per the instructions in the handout.

Computing environment

%load_ext watermark
%watermark -v -p  numpy,scipy,pandas,bokeh,colorcet,bi1x,jupyterlab

Python implementation: CPython
Python version       : 3.12.5
IPython version      : 8.27.0

numpy     : 1.26.4
scipy     : 1.13.1
pandas    : 2.2.2
bokeh     : 3.6.0
colorcet  : 3.1.0
bi1x      : 0.0.14
jupyterlab: 4.2.5