# Calculation and Interpretation of Substrate Assimilation Rates in Microbial Cells Based on Isotopic Composition Data Obtained by nanoSIMS

Lubos Polerecky, Meri Eichner, Takako Masuda, Tomáš Zavřel, Sophie Rabouille, Douglas Campbell, Kimberly Halsey

<figure><img src="image_1.png" /></figure>

## To cite this version:

Lubos Polerecky, Meri Eichner, Takako Masuda, Tomáš Zavřel, Sophie Rabouille, et al., Calculation and Interpretation of Substrate Assimilation Rates in Microbial Cells Based on Isotopic Composition Data Obtained by nanoSIMS. Frontiers in Microbiology, 2021, 12, 10.3389/fmicb.2021.621634. hal-03469894

HAL Id: hal-03469894

https://hal.science/hal-03469894

Submitted on 8 Dec 2021

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

Check for
updates

# Calculation and Interpretation of Substrate Assimilation Rates in Microbial Cells Based on Isotopic Composition Data Obtained by nanoSIMS

Lubos Polerecky¹*, Meri Eichner², Takako Masuda², Tomáš Zavřel³, Sophie Rabouille⁴,⁵, Douglas A. Campbell⁶ and Kimberly Halsey⁷*

¹ Department of Earth Sciences, Utrecht University, Utrecht, Netherlands, ² Institute of Microbiology, Czech Academy of Sciences, Centre Algatech, Treboň, Czechia, ³ Global Change Research Institute, Czech Academy of Sciences, Brno, Czechia, ⁴ Sorbonne Université, CNRS, Laboratoire d'Océanographie de Villefranche, LOV, Villefranche-sur-mer, France, ⁵ Sorbonne Université, CNRS, Laboratoire d'Océanographie Microbienne, LOMIC, Banyuls-sur-mer, France, ⁶ Department of Biology, Mount Allison University, Sackville, NB, Canada, ⁷ Department of Microbiology, Oregon State University, Corvallis, OR, United States

**OPEN ACCESS**

**Edited by:**
Qi Zhao,
University of Science and Technology
Liaoning, China

**Reviewed by:**
Celine Loussert-Fonta,
Université de Fribourg, Switzerland
Marcel Kuypers,
Max Planck Institute for Marine
Microbiology (MPG), Germany

*Correspondence:
Lubos Polerecky
l.polerecky@uu.nl
Kimberly Halsey
halseyk@science.oregonstate.edu

**Specialty section:**
This article was submitted to Systems Microbiology, a section of the journal Frontiers in Microbiology

**Received:** 27 October 2020
**Accepted:** 08 November 2021

**Published:** 30 November 2021

**Citation:**
Polerecky L, Eichner M, Masuda T, Zavřel T, Rabouille S, Campbell DA and Halsey K (2021) Calculation and Interpretation of Substrate Assimilation Rates in Microbial Cells Based on Isotopic Composition Data Obtained by nanoSIMS. *Front. Microbiol.* 12:62164. doi: 10.3389/fmicb.2021.621634

Stable isotope probing (SIP) combined with nano-scale secondary ion mass spectrometry (nanoSIMS) is a powerful approach to quantify assimilation rates of elements such as C and N into individual microbial cells. Here, we use mathematical modeling to investigate how the derived rate estimates depend on the model used to describe substrate assimilation by a cell during a SIP incubation. We show that the most commonly used model, which is based on the simplifying assumptions of linearly increasing biomass of individual cells over time and no cell division, can yield underestimated assimilation rates when compared to rates derived from a model that accounts for cell division. This difference occurs because the isotopic labeling of a dividing cell increases more rapidly over time compared to a non-dividing cell and becomes more pronounced as the labeling increases above a threshold value that depends on the cell cycle stage of the measured cell. Based on the modeling results, we present formulae for estimating assimilation rates in cells and discuss their underlying assumptions, conditions of applicability, and implications for the interpretation of intercellular variability in assimilation rates derived from nanoSIMS data, including the impacts of storage inclusion metabolism. We offer the formulae as a Matlab script to facilitate rapid data evaluation by nanoSIMS users.

**Keywords:** nanoSIMS, stable isotope probing, assimilation rates, storage inclusions, cell growth model

## INTRODUCTION

Stable isotope probing (SIP) is an experimental approach used to trace the fates of substrates within a microbial community. In this approach, samples are incubated with a target substrate that is enriched in a stable isotope, such as ¹³C, and then cells or cellular components are analyzed for enrichment in the target isotope (Boschker and Middelburg, 2002; Dumont and Murrell, 2005;
Neufeld et al., 2007; Hungate et al., 2015). SIP can be used to identify microbes in environmental samples that are actively using a particular growth substrate by selectively recovering and analyzing isotopic signatures in biomarkers such as DNA, rRNA, lipids, fatty acids or proteins (Boschker and Middelburg, 2002; Neufeld et al., 2007; Jehmlich et al., 2016). SIP can also be used to quantify rates of substrate assimilation and compartmentalization by specific community members with the goal of estimating fluxes at the ecosystem scale (Dumont and Murrell, 2005; Hungate et al., 2015). In the latter application, use of nano-scale secondary ion mass spectrometry (nanoSIMS; see Hoppe et al., 2013 or Nuñez et al., 2017 for a review of the method) to measure substrate assimilation into individual cells has become popular in the last decade (Musat et al., 2012; Pett-Ridge and Weber, 2012; Mayali, 2020; Ploug, 2021). Coupling SIP with nanoSIMS enables measurement of activities of microbes from the environment, including those that cannot yet be cultured (Krupke et al., 2015; Harding et al., 2018; Dekas et al., 2019; Geerlings et al., 2020; Mills et al., 2020), are rare (Musat et al., 2008; Zimmermann et al., 2015), or are tightly associated with their geophysical and/or biological environments (Foster et al., 2011, 2013; Milucka et al., 2012; Bonnet et al., 2016; Berthelot et al., 2019; Loussert-Fonta et al., 2020), thus offering opportunities to advance understanding of how single-celled life works, adapts, and impacts geochemical element cycling.

A typical output of a nanoSIMS measurement is the isotopic composition of a cell, determined as the cell-specific isotope ratio, *R*, or atom fraction, *x* [e.g., for carbon isotopes, *R* = $^{13}C/^{12}C$ and *x* = $^{13}C/($^{12}C+^{13}C$)]. One way of using this data is to convert *R* or *x* into quantities that can directly be interpreted in terms of the cell's mass balance, including $Fx_{net}$, $K_A$ or $X_{net}$ (Popa et al., 2007; Finzi-Hart et al., 2009; Stryhanyuk et al., 2018; Dekas et al., 2019). Specifically, the parameter $Fx_{net}$, introduced by Popa et al. (2007) and equal to the parameter $K_A$ introduced by Stryhanyuk et al. (2018), reflects the net amount of element (e.g., carbon) assimilated by the cell during the SIP experiment ($E_a$) relative to the initial content of the element in the cell prior to the SIP experiment ($E_i$), i.e., $Fx_{net} = K_A = E_a/E_i$. Similarly, the parameter $X_{net}$, used for instance by Dekas et al. (2019), reflects the net amount of element assimilated by the cell relative to the final element content in the cell ($E_f$), i.e., $X_{net} = E_a/E_f$ (Table 1).

The quantities described above do not account for incubation duration so can only be compared across SIP experiments that are incubated for the same duration. Calculating the rate of substrate assimilation obviates this requirement and has become a common output from nanoSIMS analyses of microbial cells (Ploug, 2021, and references therein). The two types of substrate assimilation rates used are the substrate-specific and cell-specific assimilation rates. While the substrate-specific rate [e.g., the rate of carbon assimilation normalized to the carbon content of the cell; in mol C (mol C)$^{-1}$ h$^{-1}$] provides useful information on cellular turnover of the substrate and can be converted to a cell's doubling time (e.g., Musat et al., 2008; Martínez-Pérez et al., 2016; Eichner et al., 2017), the cell-specific rate (e.g., the carbon assimilation rate per cell; in mol C cell$^{-1}$ h$^{-1}$) can be used to estimate the impact of a cell's population on element fluxes over larger spatial and temporal scales (e.g., Foster et al., 2011;

Krupke et al., 2015; Bonnet et al., 2016; Klawonn et al., 2016; Arandia-Gorostidi et al., 2017; Mills et al., 2020).

Calculating substrate assimilation rates by individual cells from nanoSIMS data requires consideration of cell division and the mathematical description of cell growth. Substrate assimilation during a SIP incubation, as indicated by isotopic enrichment of cell biomass at the end of the incubation (i.e., a positive value of the cell-specific $Fx_{net}$ or $X_{net}$), implies that there is a non-zero probability that the cell divided during the SIP incubation. A simple mass balance calculation reveals that, under balanced growth, this probability reaches 1 when the cell-specific quantities $Fx_{net}$ and $X_{net}$ exceed 1 and 0.5, respectively. Some recent studies have considered cell division occurring during SIP incubations by applying a correction factor (i.e., 2) in growth rate calculations (Schoffelen et al., 2018, 2019; Olofsson et al., 2019a,b). The mathematical basis for such correction factors has not been, however, explicitly formulated, and the impact of cell division on calculated rates of substrate assimilation in cells has not been investigated in a systematic way.

Calculated substrate assimilation rates will vary depending on the model used to describe substrate assimilation over the cell cycle (i.e., cell growth) and thus over the incubation time. A linear model assumes that substrate assimilation proceeds at a constant cell-specific rate (i.e., zero-order kinetics) leading to a linear increase in the biomass of individual cells over time. In contrast, an exponential model assumes that substrate assimilation proceeds at a rate that is linearly proportional to the instantaneous cell biomass (i.e., first-order kinetics) leading to an exponential increase in the biomass of individual cells over time (Collins and Richmond, 1962; Koch and Schaechter, 1962; Koch, 1966). While the linear model of substrate assimilation is most commonly used to calculate assimilation rates in cells from nanoSIMS data (e.g., Foster et al., 2013; Klawonn et al., 2016; Stryhanyuk et al., 2018; Mills et al., 2020), the exponential model is used rather rarely (e.g., Martínez-Pérez et al., 2016; Berthelot et al., 2019; Geerlings et al., 2020, 2021; Polerecky et al., 2021). The impact of the substrate assimilation model (linear vs. exponential) on the calculated assimilation rate has previously not been studied in detail.

Assimilation rates determined at single-cell resolution often reveal intercellular heterogeneity in populations. Previous nanoSIMS studies typically attributed such heterogeneity to factors that affect the *assimilation* of an (isotopically labeled) element during the incubation (causing variation in $E_a$), i.e., differences in the intrinsic metabolic activity of a cell, which may vary with cell cycle stage or cell age, or stochastic gene expression (e.g., reviewed by Ackermann, 2015). Other sources of heterogeneity have received less attention, such as assimilation of unlabeled sources of the target element, which will affect the degree of isotopic labeling of the assimilated material, or differences in the cell's initial elemental content (variation in $E_i$). For instance, selective assimilation of C and N into storage inclusions or mobilization of C or N from existing cell material during a SIP experiment can critically affect both the (apparent) amount of an element assimilated during an incubation ($E_a$) and the initial amount of an element present in the cell ($E_i$) (Polerecky et al., 2021). Thus, metabolism of storage compounds can partly
<table><caption>**TABLE 1** | List of symbols used in this study.</caption><thead><tr><th>Symbol</th><th>Unit</th><th>Description</th></tr></thead><tbody><tr><td>C, C<sub>i</sub></td><td>mol C cell<sup>-1</sup></td><td>Carbon content of a cell; 'i' refers to the initiation of the SIP incubation, i.e., time-point t = 0.</td></tr><tr><td>C<sub>max</sub></td><td>mol C cell<sup>-1</sup></td><td>Carbon content of a cell just before binary cell division; C<sub>max</sub>/2 corresponds to the carbon content just after binary cell division.</td></tr><tr><td>s</td><td>-</td><td>Cell cycle stage, calculated from the C content as s = C/(C<sub>max</sub>/2) - 1.</td></tr><tr><td>⟨C⟩</td><td>mol C cell<sup>-1</sup></td><td>Carbon content of a cell averaged over the cell cycle. It is equal to the average in a population of cells with perfectly unsynchronized cell cycles. It is related to C<sub>max</sub> according to ⟨C⟩ = C<sub>max</sub>/(2·ln(2)) and ⟨C⟩ = C<sub>max</sub>·ln(2) for cells assimilating C with zero-order and first-order kinetics, respectively (Koch, 1966).</td></tr><tr><td>x, x<sub>i</sub></td><td>-</td><td><sup>13</sup>C atom fraction of a cell, defined as x = <sup>13</sup>C/(<sup>13</sup>C+<sup>13</sup>C); 'i' refers to t = 0.</td></tr><tr><td>x<sub>S</sub><sup>E</sup></td><td>-</td><td>Source-normalized excess <sup>13</sup>C atom fraction of a cell (Eq. 10).</td></tr><tr><td>r</td><td>mol C cell<sup>-1</sup> h<sup>-1</sup></td><td>Cell-specific rate of carbon assimilation by a cell.</td></tr><tr><td>⟨r⟩</td><td>mol C cell<sup>-1</sup> h<sup>-1</sup></td><td>Average cell-specific C assimilation rate in a population.</td></tr><tr><td>k</td><td>h<sup>-1</sup></td><td>Carbon-specific rate of carbon assimilation by a cell.</td></tr><tr><td>t</td><td>h</td><td>Doubling time; time needed for a cell to double its carbon content. Calculated as τ = C<sub>max</sub>/(2·r) and τ = (ln2)/k for a cell assimilating C with zero-order and first-order kinetics, respectively.</td></tr><tr><td>ρ</td><td>mol C µm<sup>-3</sup></td><td>Carbon density of a cell.</td></tr><tr><td>x<sub>S,tar</sub></td><td>-</td><td><sup>13</sup>C atom fraction of the target source provided externally during the SIP incubation.</td></tr><tr><td>x<sub>S,alt</sub></td><td>-</td><td><sup>13</sup>C atom fraction of an alternative carbon source; can be external (e.g., present in the environment of the cell) or internal (present in the form of intra-cellular C storage inclusions).</td></tr><tr><td>x<sub>S,eff</sub></td><td>-</td><td>Effective <sup>13</sup>C atom fraction of the carbon source (Eq. 6).</td></tr><tr><td>f<sub>tar</sub></td><td>-</td><td>Fraction of C assimilated by the cell from the target source; the remaining fraction of C is assimilated from the alternative C source (f<sub>alt</sub> = 1 - f<sub>tar</sub>).</td></tr><tr><td>θ</td><td>-</td><td>Probability density function (PDF) describing the distribution of C among cells in a population. Given by Eq. 16 and 18 for a population with perfectly unsynchronized and partially synchronized cell cycles, respectively. Examples shown in **Figure 4B**, **Supplementary Figures 3 and 4**.</td></tr><tr><td>ζ</td><td>-</td><td>Probability density function (PDF) describing the distribution of x<sub>S</sub><sup>E</sup> in a population of cells assimilating C with zero-order kinetics. Examples shown in **Figure 4C**.</td></tr><tr><td>γ</td><td>-</td><td>Degree of cell cycle synchronicity characterizing the distribution of C content among cells in a population with partially synchronized cell cycles. Examples shown in **Supplementary Figures 4A,B**.</td></tr><tr><td>X<sub>net</sub></td><td>-</td><td>Net amount of element assimilated by the cell (E<sub>a</sub>) relative to the final element content in the cell (E<sub>i</sub>). Calculated according to X<sub>net</sub> = E<sub>a</sub>/E<sub>i</sub>. Used, for instance, by Dekas et al. (2019). For a cell that did not divide during the SIP incubation, X<sub>net</sub> = x<sub>S</sub><sup>E</sup> (Supplementary Material, Section "Relationships Between x<sub>S</sub><sup>E</sup>, X<sub>net</sub>, FX<sub>net</sub> and K<sub>A</sub>").</td></tr><tr><td>FX<sub>net</sub>, K<sub>A</sub></td><td>-</td><td>Net amount of element assimilated by the cell (E<sub>a</sub>) relative to the initial element content in the cell (E<sub>i</sub>). Calculated according to FX<sub>net</sub> = K<sub>A</sub> = E<sub>a</sub>/E<sub>i</sub>. Parameter FX<sub>net</sub> was introduced by Popa et al. (2007), while parameter K<sub>A</sub> was introduced by Stryhanyuk et al. (2018). For a cell that did not divide during the SIP incubation, FX<sub>net</sub> = K<sub>A</sub> = x<sub>S</sub><sup>E</sup>/(1 - x<sub>S</sub><sup>E</sup>) (Supplementary Material, Section "Relationships Between x<sub>S</sub><sup>E</sup>, X<sub>net</sub>, FX<sub>net</sub> and K<sub>A</sub>").</td></tr></tbody></table>

explain variability in the enrichment of cellular biomass (e.g., FX<sub>net</sub> or X<sub>net</sub>).

In this work we use mathematical modeling to identify key processes and parameters that play a role when calculating and interpreting cellular substrate assimilation rates from the isotopic composition information determined by SIP-nanoSIMS. First, we investigate how the isotopic composition of a cell will vary in time depending on the model of substrate assimilation (zero-order vs. first-order kinetics), the substrate assimilation rate, cell biomass, and the isotopic composition of the substrate. We specifically analyze how cell division during the incubation impacts the predicted isotopic composition of a cell and how the isotopic compositions will vary among cells within a population depending on the degree of cell cycle synchronicity. Based on these theoretical analyses, we provide nanoSIMS users with step-by-step guidelines for calculating substrate assimilation rates including formulas and their underlying assumptions. We also identify key biological factors that can introduce variability into the calculated rates, specifically highlighting the previously unrecognized role of storage inclusions and their impact on interpretations of intercellular variability in calculated assimilation rates.

## MATERIALS AND METHODS

Our analysis is based on mass balances, which we formulate as differential equations to facilitate identification of parameters impacting the predicted isotopic composition of cells during a SIP experiment. We focus our analysis on assimilation of carbon (C) from a <sup>13</sup>C-labeled carbon source. Analogous analyses can be conducted for other elements and their stable isotopes commonly used in SIP experiments (e.g., <sup>15</sup>N, <sup>18</sup>O). Throughout this study, we generally adhere to the accepted notation guidelines for reporting isotope enrichment data (Coplen, 2011) with a few exceptions (see **Table 1** for a list of symbols and their definition used in this study).

### <sup>13</sup>C Labeling of an Individual Cell

The model of C assimilation by an individual cell is mathematically formulated by a differential equation

$$ \frac{dC(t)}{dt} = r(t), \qquad (1) $$

where C denotes the cell-specific C content (in mol C cell<sup>-1</sup>) and r denotes the cell-specific C assimilation rate (in mol C
cell⁻¹ h⁻¹). Both quantities are a function of time, as indicated by the argument (t). In this study, we neglect the potential effect of kinetic isotopic fractionation as it is typically minor compared to the changes in the isotopic composition of the cell due to labeled substrate assimilation during the SIP incubation. If needed, the effect can easily be included in the analysis described below. Additionally, we do not consider processes associated with cellular biomass turnover, e.g., exudation of freshly fixed C or respiration of storage compounds. Such processes could be implemented by assuming a second, “loss” term in Eq. 1. Their mathematical treatment is, however, beyond the scope of this study.

Given the rate of increase in the total C (i.e., ¹²C + ¹³C) content of the cell, r, the ¹³C content of the cell will increase at a rate r multiplied by the probability that the assimilated C atom is ¹³C. This probability is equal to the ¹³C atom fraction of the carbon source, denoted as xS. Consequently, the ¹³C content of the cell is described by a differential equation

$$ \frac{d^{13}C(t)}{dt} = r(t) \cdot x_S. \qquad (2) $$

To make our analysis more general, we assume that the assimilated C may originate from two sources: (1) the target C source, which is isotopically labeled by the experimenter by the addition of a naturally less abundant stable isotope (e.g., ¹³C), and (2) an alternative C source, which only contains the natural abundance of ¹³C and may be present and assimilated by the cell without experimenter's knowledge. For example, ¹³C-labeled dissolved inorganic carbon (DIC) is a common target C source and unlabeled dissolved organic carbon (DOC) could be an alternative C source. Assuming that C assimilations from the target and alternative sources are independent of each other, the total C and ¹³C contents of the cell will vary with time according to differential equations

$$ \frac{dC(t)}{dt} = r_{tar}(t) + r_{alt}(t), \qquad (3) $$

$$ \frac{d^{13}C(t)}{dt} = r_{tar}(t) \cdot x_{S,tar} + r_{alt}(t) \cdot x_{S,alt}, \qquad (4) $$

where subscripts “tar” and “alt” refer to the target and alternative C source, respectively.

Although in some SIP experiments the isotopic composition of the C sources may vary in time (e.g., Schreiber et al., 2016; Geerlings et al., 2020), we consider hereafter the more common situation that xS,tar and xS,alt are time-independent. We make this assumption because our aim is to derive simple formulas to calculate the substrate assimilation rate based on the isotopic composition of the cell. A more advanced analysis describing SIP experiments with time-dependent substrate labeling is provided in Supplementary Methods, Section “Modeling Cellular Assimilation of Substrates With Time-Dependent Isotopic Composition,” and in Supplementary Figures 1,2. Considering this assumption, we denote the total cell-specific C assimilation rate as r = rtar + ralt, and the fraction of C assimilated by the cell from the target source as

fstar = rstar/(rstar + ralt), which makes Eq. 3 identical to Eq. 1 and yields the following differential equation for the cellular ¹³C content:

$$ \frac{d^{13}C(t)}{dt} = r(t) \cdot x_{S,eff}, \qquad (5) $$

where

$$ x_{S,eff} = f_{tar} \cdot x_{S,tar} + (1 - f_{tar}) \cdot x_{S,alt} \qquad (6) $$

is the effective ¹³C atom fraction of the carbon source.

To understand the consequences of C assimilation from a source characterized by xS,eff for the isotopic composition of the cell requires the solutions of differential equations 1 and 5. The result depends on how r may or may not vary across the cell cycle. We consider two idealized scenarios: (1) r is constant across the cell cycle (i.e., the cell assimilates C according to zero-order kinetics), (2) r is linearly proportional to the instantaneous cell-specific C content, r = k · C (i.e., the cell assimilates C according to first-order kinetics), where k is constant across the cell cycle. The model of zero-order kinetics can be conceptualized as a cell whose C content increases due to the constant activity of the cell components (e.g., ribosomes, proteins) that are present at the beginning of the cell cycle and do not change in concentration over its course (Figure 1A). The model of first-order kinetics is conceptualized as described by Koch and Schaechter (1962) where the cell components are synthesized continuously and operate at constant rates throughout the cell cycle such that the ratio between these operational units and the mass of the total cell remains constant across the cell cycle (Figure 1B). Additionally, we assume that the C density of a cell, denoted hereafter as ρ (in mol C µm⁻³), is constant across the cell cycle. Variation in ρ will occur if C is assimilated into storage inclusions (e.g., polysaccharide or lipid inclusions) or if C storage inclusions

C assimilation by a cell: dC/dt = r

A Zero-order kinetics: r = constant

<figure><img src="image_4.png" /></figure>

B First-order kinetics: r = k·C, k = constant

<figure><img src="image_5.png" /><figcaption>**FIGURE 1** | Conceptual diagram of C assimilation during cell growth. Cell is represented by a rectangle with round edges. Size of the rectangle represents C content of the cell, thickness of the arrow entering the cell represents the cell-specific rate of C assimilation, r. (A) Zero-order kinetics of C assimilation assumes that r is independent of the cellular C content. (B) First-order kinetics of C assimilation assumes that r increases proportionally to the instantaneous cellular C content.</figcaption></figure>
are catabolized for respiration or the synthesis of other cell components during certain periods of the cell cycle. Effects of these variations are qualitatively examined in the Discussion.

## Zero-Order Kinetics of C Assimilation

Assuming constant r, solutions to the differential equations 1 and 5 are readily found as

$$C(t) = C_i + r \cdot t, \quad (7)$$

$$^{13}C(t) = C_i \cdot x_i + r \cdot x_{S,eff} \cdot t, \quad (8)$$

where Cᵢ and xᵢ denote, respectively, the C content and ¹³C atom fraction of the cell at the time-point when the SIP experiment with the labeled substrate was initiated (i.e., at t = 0). Thus, zero-order kinetics of C assimilation corresponds to a linear increase in cell biomass.

Equations 7-8 imply that the ¹³C atom fraction of the cell, defined as x = ¹³C/(¹²C+¹³C), increases in time according to

$$x(t) = \frac{^{13}C(t)}{^{12}C(t) + ^{13}C(t)} = \frac{^{13}C(t)}{C(t)} = \frac{C_i \cdot x_i + r \cdot x_{S,eff} \cdot t}{C_i + r \cdot t}. \quad (9)$$

To simplify this function, we define the source-normalized excess ¹³C atom fraction as

$$x_S^E(t) = \frac{x(t) - x_i}{x_{S,eff} - x_i}. \quad (10)$$

This is a convenient quantity for describing the isotopic composition of the cell because it increases from 0 to 1 as the cell increases its ¹³C atom fraction from xᵢ towards x_{S,eff} (see Supplementary Methods, Section “Relationships between x_S^E, X_{net}, Fx_{net} and K_A,” for a relationship between x_S^E and the quantities Fx_{net}, K_A and X_{net} introduced previously by other authors). Upon rearrangement, this definition combined with Eq. 9 yields

$$x_S^E(t) = \frac{r \cdot t}{C_i + r \cdot t}. \quad (11)$$

The next step is to account for cell division. Specifically, we assume binary division, in which the cell divides into two identical cells when its C content reaches some critical value C_{max} (Collins and Richmond, 1962; Koch, 1966; Westfall and Levin, 2017). Hence, an important assumption underlying the models in this study is that cells are not starving, a condition that is commonly associated with accumulation of storage products with no regular cell division. Thus, equations 7-11 correctly describe the C content and ¹³C atom fraction of a cell if, and only if, C(t) is in the interval between C_{max}/2 and C_{max}. When C(t) in equations 7-11 reaches C_{max}, it must be reset to C_{max}/2. As a result, ¹³C labeling of the cell, described by the parameter x_S^E, will vary in time following a “zig-zag” pattern that depends on Cᵢ, C_{max}, and r (see Results, Figure 2A). Hereafter, we denote the function describing this pattern as Z(t, Cᵢ, C_{max}, r). The consequences of the zig-zag pattern on the resultant substrate assimilation rate are discussed in Results. Note that although zero-order kinetics of C assimilation implies a linear increase in the biomass of an individual cell (Eq. 7), binary cell division of such cells yields a

population growing exponentially in cell number and total cell biomass over time.

## First-Order Kinetics of C Assimilation

Under first-order kinetics, r is linearly proportional to the instantaneous cell-specific C content, i.e., r(t) = k · C(t), where k is a constant. The parameter k is commonly referred to as the growth rate constant (in h⁻¹) and represents the carbon-specific rate of C assimilation [in mol C (mol C)⁻¹ h⁻¹]. Substituting this rate expression into Eq. 1, we obtain for the cell-specific C content

$$C(t) = C_i \cdot e^{k \cdot t}. \quad (12)$$

where Cᵢ is the initial C content of the cell. Thus, first-order kinetics of C assimilation corresponds to an exponentially growing cell. Note that in some studies (e.g., Martínez-Pérez et al., 2016) exponential cell growth is assumed to be described by a function 2^Growthrate·t rather than e^{k·t}. These descriptions are equivalent when “Growthrate” is calculated as k/ln(2). The corresponding doubling time, τ, is calculated by considering that C(τ) = 2·Cᵢ, which yields a formula τ = ln(2)/k = 1/Growthrate.

To evaluate the dynamic of the ¹³C atom fraction (x) for an exponentially growing cell, we start from the definition of x and apply the quotient rule to evaluate dx/dt. Then, we use equations 1, 5 and 12 yielding the following differential equation for x:

$$\frac{dx}{dt} = -k \cdot [x - x_{S,eff}]. \quad (13)$$

The solution to this differential equation is

$$x(t) = x_{S,eff} - [x_{S,eff} - x_i] \cdot e^{-k \cdot t}. \quad (14)$$

Using the definition in Eq. 10, this function is simplified to

$$x_S^E(t) = 1 - e^{-k \cdot t}. \quad (15)$$

when the isotopic composition of the cell is described by x_S^E instead of x. Thus, ¹³C labeling of a cell assimilating C according to first-order kinetics is described by one minus an exponential function of time.

To account for cell division, we again assume that the cell divides into two identical cells when its C content reaches some critical value C_{max}. Thus, when C(t) described by Eq. 12 reaches C_{max}, it must be reset to C_{max}/2. In contrast to cell-specific C assimilation modeled by zero-order kinetics, this reset in C(t) due to cell division does not affect the time-dependence of the isotopic composition of the cell assimilating C according to first-order kinetics. This is because the assimilation rate decreases by half as the cell divides, hence the ratio of ¹³C uptake to ¹²C present in the cell remains the same. This insight also emerges directly from the differential equation 13, which shows that x, and hence also x_S^E, is independent of Cᵢ (see Eq. 14-15).

## ¹³C Labeling in a Population of Cells

Here, we expand our analysis from the single cell view to a population of cells with varying degrees of cell cycle synchronization. Our aim is to reveal how the predicted isotopic
<figure><img src="image_6.png" /><figcaption>**FIGURE 2** | Simulations of the time-dependence of the isotopic composition and C content in cells assimilating C according to different models. Shown are the source-normalized excess $^{13}$C atom fraction (panels A,D), total C content (panels B,E), and $^{13}$C content (panels C,F) for two cells with different initial C contents ($C_i = C_{max}/2$ for cell 1, $C_i = 1.5 \cdot C_{max}/2$ for cell 2). Note that cells divide when their total C content reaches $C_{max}$. For comparison, results calculated based on a model that ignores cell division are shown by dashed lines in panel (A). The number of primes behind the digit 1 and 2 indicates the generation of the daughter cells corresponding to the original mother cells. The cells assimilating C with zero-order kinetics had equal cell-specific C assimilation rates, r, while the cells assimilating C with first-order kinetics had equal carbon-specific C assimilation rates, k. These rates were related as $k = r/\langle C \rangle$, where \langle C \rangle was the same for both cell types. In all panels the time was normalized by the doubling time $\tau$, calculated as $\tau = C_{max}/(2 \cdot r)$ and $\tau = \ln(2)/k$ for the cells assimilating C with zero-order and first-order kinetics, respectively. Blue arrows with numbers at $t = n \cdot \tau$ indicate $x_S^E = 1 - (1/2)^n$, where n = 1, 2, etc. Simulations were performed for $x_i = 0.011$ and for hypothetical values of $x_{S,eff} = 0.2$, $\langle C \rangle = 10$ fmol C cell$^{-1}$ and $r = 1.35$ fmol C cell$^{-1}$ h$^{-1}$. These values correspond to $k = 0.135$ h$^{-1}$, $\tau \approx 5.1$ h, and $C_{max} \approx 13.9$ and $14.4$ fmol C cell$^{-1}$ for the cells assimilating C with zero-order and first-order kinetics, respectively.</figcaption></figure>

enrichment varies within a population of cells assimilating C at equal rates and how this variability depends on the model used to describe C assimilation by the cell and the synchronization of the cell cycles among cells. Note that regardless of whether C assimilation is modeled by zero- or first-order kinetics, binary division at the end of each cell cycle causes the population to grow exponentially. Moreover, the population doubling time, $\tau$, is related to the C assimilation rate according to $\tau = C_{max}/(2 \cdot r)$ and $\tau = (\ln 2)/k$ for the models assuming zero-order and first-order kinetics of C assimilation, respectively.

## Zero-Order Kinetics of C Assimilation in a Population

First, we consider zero-order kinetics of C assimilation by individual cells. If cell cycles in a population are not synchronized, C at a given point in time varies among cells between $C_{max}/2$ and $C_{max}$ (Collins and Richmond, 1962; Koch, 1966). Thus, because of the dependence of the Z function on $C_i$ (see above, Section “Zero-order Kinetics of C Assimilation”), individual cells will have different $x_S^E$ depending on their $C_i$ even if they assimilate C at the same rate r. We denote the probability density function (PDF) describing the distribution of $C_i$ and $x_S^E$ among cells by $\theta$ and $\zeta$, respectively, and first consider perfectly

unsynchronized and then partially synchronized cells. In both cases we assume that the critical C content where each cell divides, $C_{max}$, is equal for all cells.

As shown by Koch (1966), the PDF for perfectly unsynchronized cells is

$$ \theta(m) = \frac{8 \cdot \ln(2)}{C_{max}} \cdot e^{-\frac{2 \cdot m \cdot \ln(2)}{C_{max}}} \quad (16) $$

for m in the interval $C_{max}/2 \le m \le C_{max}$, and it is equal to zero for m outside this interval. For these cells, the C content of the average cell in the population, which is equivalent to the C content averaged over the cell cycle, is related to the critical C content $C_{max}$ by the following formula (Koch, 1966):

$$ \langle C \rangle = \int_{0}^{\infty} m \cdot \theta(m) dm = \frac{C_{max}}{2 \cdot \ln(2)} \approx 0.72 \cdot C_{max}. \quad (17) $$

Furthermore, both $\theta(m)$ and $\langle C \rangle$ are time independent as the population grows.

For partially synchronized cells, $\theta$ can take various forms. Here, we assume $\theta$ to be based on a Gaussian function. However, because $\theta$ is non-zero within the interval between $C_{max}/2$ and $C_{max}$, the tails of the Gaussian function reaching outside this
interval must be constrained within this interval with a factor that accounts for binary cell division (Supplementary Figure 3). This consideration yields

$$ \theta(m) = A \cdot \sum_{j = -\infty}^{\infty} 2^j \cdot e^{-\frac{(m + j \cdot C_{max}/2 - C_0)^2}{2(\Delta C)^2}} \quad (18) $$

for $m$ in the interval $C_{max}/2 \le m \le C_{max}$, and zero for $m$ outside this interval. In this function, $C_0$ and $\Delta C$ describes the center and width of the Gaussian function, respectively ($C_0$ must lie between $C_{max}/2$ and $C_{max}$), and $A$ is the normalization constant such that $\int_{C_{max}/2}^{C_{max}} \theta(m) dm = 1$.

We define the degree of cell cycle synchronicity in the population by the ratio $\gamma = C_{max}/(8 \cdot \Delta C)$ (Supplementary Figure 3). Thus, a greater degree of synchronicity corresponds to a smaller $\Delta C$, i.e., a narrower distribution of $C$ among cells (Supplementary Figures 4A,B), and vice versa. Specifically, at the limits of $\Delta C \to 0$ or $\Delta C \to \infty$, the function $\theta(m)$ in Eq. 18 describes populations of cells with perfectly synchronized ($\gamma \to \infty$) or perfectly unsynchronized ($\gamma \to 0$) cell cycles, respectively. In the latter case, $\theta(m)$ in Eq. 18 and Eq. 16 become equivalent as required (Supplementary Figure 4C). Note that if $\gamma$ is large (i.e., $\Delta C \ll C_{max}$) and each cell in the population grows at the same cell-specific rate $r$, Eq. 18 implies that 95% of cells divide within a time interval $\Delta t_{95} = 4 \cdot \Delta C/r$ (Supplementary Figure 3). Thus, for a population with cell cycles that are highly synchronized, the parameter $\gamma$ is equal to the ratio between the population doubling time, $\tau = C_{max}/(2 \cdot r)$ (Table 1), and the interval during which 95% of cells undergo division (i.e., $\gamma = \tau/\Delta t_{95}$; Supplementary Figure 3).

We performed Monte-Carlo simulations in Matlab to evaluate $\zeta$ from $\theta$. Specifically, a cell $j$ with a random initial C content $C_{ij}$ is selected based on the PDF $\theta$. Then, for a given value of $r$ and $t$, the isotopic composition of daughter cells originating from cell $j$ is calculated as $x_{Sj}^E = Z(t, C_{ij}, C_{max}, r)$. Additionally, the number of daughter cells originating from cell $j$ is calculated from the number of cell divisions as $N_j = 2^{\lfloor n_j \rfloor}$, where $\lfloor n_j \rfloor$ denotes the largest integer that is smaller than $n_j = (C_{ij} + r \cdot t)/(C_{max}/2) - 1$. Based on the values of $x_{Sj}^E$ and $N_j$ obtained for many random choices of $C_{ij}$, a histogram approximating the PDF $\zeta$ is reconstructed (see Results, Figure 4C). Finally, the average $x_S^E$ for all cells in the population is calculated as

$$ \langle x_S^E \rangle = \sum_k x_{S,k}^E \cdot \zeta(x_{S,k}^E), \quad (19) $$

where the index $k$ refers to the $k^{th}$ bin in the histogram approximating $\zeta$. Note that both the function $\zeta$ and the average value $\langle x_S^E \rangle$ depend on time, $t$, although this dependence was omitted in this expression to simplify notation.

## First-Order Kinetics of C Assimilation in a Population

Equations 14-15 show that isotopic composition of a cell assimilating C according to first-order kinetics is independent of $C_i$. Thus, in contrast to a cell assimilating C according to zero-order kinetics, the average $x_S^E$ for the population of cells

assimilating C by first-order kinetics is described by the same equation as for an individual cell (Eq. 15), and the variance among cells is zero.

## RESULTS

### Modeling ¹³C Labeling in Cells

Equations 9-11 imply that ¹³C labeling of a cell assimilating C according to zero-order kinetics depends on the cell-specific rate of C assimilation, $r$, and the initial C content of the cell, $C_i$. We illustrate the consequences of these dependencies by considering two cells with identical $r$ but different $C_i$. We assume $C_{i1} = C_{max}/2$ for cell 1, which corresponds to a cell immediately after division (initial cell cycle stage $s_{i1} = 0$), and $C_{i2} = 1.5 \cdot C_{max}/2$ for cell 2, which corresponds to a cell in the middle of its cell cycle ($s_{i2} = 0.5$; see Table 1 for the definition of $s$). Both cells then assimilate C at the same and constant rate $r$ over three cell cycles (Figures 2A-C). The lower $C_{i1}$ compared to $C_{i2}$ causes $x_{S1}^E$ to initially increase faster than $x_{S2}^E$ (Figure 2A, $t/\tau = \ll 0.5$) and the C content of cell 2 to reach $C_{max}$ earlier than that of cell 1. At $C_{max}$, the rate of increase in $x_S^E$ in the two daughter cells of cell 2 abruptly increases 2-fold due to the abrupt decrease in their C content from $C_{max}$ to $C_{max}/2$ and an unchanged $r$ (Figures 2A,B, $t/\tau = 0.5$). A similar “kink” in the evolution of $x_S^E$ is observed after cell 1 divides, but it occurs later because of the lower initial C content of cell 1 (Figure 2A, $t/\tau = 1$).

The $x_S^E$ values for the two cells described in Figure 2 become equal at regular time intervals. These intervals are separated by $\tau = C_{max}/(2 \cdot r)$, i.e., the time needed for each cell to double its C content (Figure 2A). This result is expected because the total C content of the daughter cells at time-points $t = n \cdot \tau$ ($n = 1, 2, \dots$) will be the same as the total C content of the original mother cell at $t = 0$ (Figures 2B,C). Hence, the zero-order kinetics model of C assimilation yields $x_S^E = 1 - (1/2)^n$ at time-points $t = n \cdot \tau$ irrespective of the value of $C_i$ (Figure 2A).

The time-dependence of $x_S^E$ in a cell assimilating C according to first-order kinetics is described by a function that is independent of $C_i$ (Eq. 15). Setting the initial C contents of two cells to $C_{max}/2$ and $1.5 \cdot C_{max}/2$ leads to different time-dependences of the total C and ¹³C contents for each cell (Figures 2E,F). However, the increase in $x_S^E$ over time is the same for both cells (Figure 2D). Thus, cell division and the cell's initial C content do not influence how the cell's isotopic composition varies in time if C assimilation can be adequately captured by the first-order kinetics model.

While $x_S^E$ of an individual cell modeled by zero-order kinetics is described by a function with abruptly changing time-derivative (Figures 2A, 3A), $x_S^E$ for the average cell representing its population is described by a smooth function. If cell cycles in the population are perfectly unsynchronized, the function is identical to the function describing $x_S^E$ of a cell modeled by first-order kinetics with $k = r/\langle C \rangle$, where $\langle C \rangle$ is the average C content of the cells across the cell cycle (compare solid and dashed lines in Figure 3B). Thus, for an average cell in a population with perfectly unsynchronized cell cycles, both the C
<figure><img src="image_7.png" /><figcaption>**FIGURE 3** | Monte-Carlo simulations of the time-dependence of $x_S^E$ in a population of cells with perfectly unsynchronized cell cycles. All cells assimilate C at the same cell-specific rate, $r$. (A) Time-dependence of $x_S^E$ in cells randomly selected from the population (different colors correspond to cells with different initial C content). $x_S^E$ for a dividing cell follows a pattern described by the zig-zag function Z (see **Figure 2**). For comparison, results calculated based on a model that ignores cell division are shown by dashed lines. (B) Time-dependence of $x_S^E$ for the average cell assimilating C with zero-order kinetics (black solid line), and for a cell that assimilates C with first-order kinetics with the carbon-specific rate $k = r/\langle C \rangle$ (red dashed line overlapping the black solid line). Gray dash-dotted lines show the minimum and maximum values of $x_S^E$ for cells assimilating C with zero-order kinetics. The corresponding width of the interval between the minimum and maximum values of $x_S^E$ is shown in panel (C). Simulations were performed for the same values of $\langle C \rangle$ and $r$ as those shown in **Figure 2**.</figcaption></figure>

content and ¹³C labeling behave in the same way regardless of whether their C assimilation is described by zero-order or first-order kinetics. This important result underpins the procedure for calculating rates of C assimilation in cells based on nanoSIMS data (see next section).

If cell cycles in a population are partially synchronized (as is typical, e.g., for photoautotrophs that tune cell division according to day-night cycles), the average $x_S^E$ predicted by the zero-order kinetics model follows a zig-zag pattern like that describing $x_S^E$ for a single cell (Figure 4A). This behavior is a result of the variation of the probability density function $\zeta$ becoming increasingly pronounced with the increasing degree of cell cycle synchronicity (Figure 4C). However, the interval of $x_S^E$ where $\zeta$ is non-zero is independent of the degree of cell cycle synchronicity (Figure 4C). Thus, the value of $x_S^E$ for any cell from a population with partially synchronized cell cycles, including the cell representing the average of the population, always lies within the interval defined by the minimum and maximum value of $x_S^E$ determined for cells with perfectly unsynchronized cell cycles (Figures 4A, 4C). Note that the width of this interval varies in time following a wax-and-wane pattern oscillating between zero and a maximum value that progressively decreases with time (Figure 3C). This pattern is a consequence of binary division and the assumption that the critical C content when cells divide, $C_{max}$, is equal among cells.

Together, these results show that cell division causes the time-dependence of $x_S^E$ in a cell modeled by zero-order kinetics to follow a zig-zag pattern that closely follows the exponential function in Eq. 15 describing the time-dependence of $x_S^E$ in a cell modeled by first-order kinetics (Figures 2A, 3A,B, 4A). In contrast, the assumption of zero-order kinetics of C assimilation with no accounting for cell division leads to a dramatically different pattern, where $x_S^E$ approaches 1 at a much slower rate and the variability among cells, caused for example by intercellular heterogeneity in the C content at the time of the labeled substrate addition, remains pronounced over longer time scales (compare dashed and solid lines in Figure 3A).

## Calculating Rates of Cellular C Assimilation

Using the findings presented above, we suggest a three-step procedure for calculating the cellular rates of C assimilation based on ¹³C atom fractions measured by nanoSIMS (Figure 5). We emphasize the underlying assumptions of each step to provide the foundations of the formulas. Numerical implementation of the calculation steps, including the calculation of the best estimate and uncertainty of the rates, is provided as a Matlab script that can be interfaced with experimental data organized in a spreadsheet. The script is available at https://github.com/lpolerecky/LARS. Previous studies documented that isotopic composition of cells can be significantly affected by label loss or dilution associated with common techniques of sample preparation for nanoSIMS analysis including chemical fixation, dehydration and resin embedding (Musat et al., 2012; Hoppe et al., 2013; Loussert-Fonta et al., 2020). If needed, these effects can be corrected for as previously described (Stryhanyuk et al., 2018) before applying the steps described below to calculate assimilation rates.

Step 1 converts the measured ¹³C atom fraction of the cell, $x$, to the source-normalized excess ¹³C atom fraction, $x_S^E$ (Eq. 10). This step accounts for variability in ¹³C labeling of the cell due to the ¹³C-labeled source, $x_{S,eff}$ (see Eq. 6). Quality results require that $x_{S,eff}$ be well constrained by the experimental set-up, such as by direct measurement. The initial ¹³C atom fraction of the cell, $x_i$, is obtained by nanoSIMS measurements of control (unlabeled) cells that were prepared in the same manner as treated (labeled) cells.

Step 2 calculates the carbon-specific rate of C assimilation, $k$, according to

$$k = -\frac{1}{t} \cdot \ln(1-x_S^E), \qquad (20)$$

where $x_S^E$ is obtained in Step 1. This calculation accounts for variability in ¹³C labeling of the cell due to the incubation time (Eq. 15). The underlying assumption of this step is that the measured cell assimilated C during the SIP incubation according
<figure><img src="image_8.png" /><figcaption>**FIGURE 4** | Monte-Carlo simulations of the time-dependence of $x_S^E$ in a population of cells with partially synchronized cell cycles. All cells assimilated C at the same cell-specific rate, *r*. Simulations were performed for the same values of $\langle C \rangle$ and *r* as those shown in **Figure 2**. The distribution of $C_i$ among cells is described by Eq. 18 with $C_0 = 9$ fmol C cell$^{-1}$ (red) or $C_0 = 13$ fmol C cell$^{-1}$ (blue) and the degree of synchronicity $\gamma = 1.73$ (corresponding to $\Delta C = 1$ fmol C cell$^{-1}$, both cases). Results for perfectly unsynchronized cells ($\gamma = 0$) are shown for comparison (black). (A) Time-dependence of the average $x_S^E$. Dash-dotted lines show the minimum and maximum values of $x_S^E$ for a population of perfectly unsynchronized cells (reproduced from **Figure 3B**). (B) PDF describing the distribution of $C_i$ among cells. (C) PDF describing the distribution of $x_S^E$ among cells, shown at two time-points ($0.4 \tau$ and $1.6 \tau$, where $\tau$ is the population doubling time; see legend). Colors and symbols in panels (A–C) correspond to each other.</figcaption></figure>

to first-order kinetics or, equivalently, that the measured cell represents the average cell in a population that assimilated C according to zero-order kinetics, had the same cell-specific rate *r*, and had perfectly unsynchronized cell cycles (Figure 3B). Importantly, this step accounts for the possibility that the measured cell is the product of cell division that occurred during the SIP incubation.

Step 3 calculates the cell-specific rate of C assimilation, *r*. Because *r* is an absolute measure of the C assimilation rate (mol C cell$^{-1}$ h$^{-1}$), its accuracy depends on the knowledge of the C content of the measured cell (Supplementary Methods, Section “Estimating Cellular C Content”). There are three approaches to calculate *r*. Approach A applies when the C content of the measured cell is not known precisely but can be approximated by the average C content of the cell over a cell cycle, $\langle C \rangle$ (Figure 5, Step 3A). In this case, *r* is calculated according to

$$r = k \cdot \langle C \rangle \qquad (21)$$

using the value of *k* obtained in Step 2. Approach B applies when only the C content of the measured cell, *C*, but not the average C content $\langle C \rangle$, can be constrained (Figure 5, Step 3B). In this case, *r* is calculated according to

$$r = k \cdot C \qquad (22)$$

using the value of *k* obtained in Step 2. Finally, Approach C applies when both *C* and $\langle C \rangle$ can be constrained (Figure 5, Step 3C). In this case, *r* is calculated according to

$$r = Z^{-1}(t, C, C_{max}, x_S^E) \qquad (23)$$

using the value of $x_S^E$ obtained in Step 1. Here, $Z^{-1}$ is the inverse of the zig-zag function introduced above (Figure 2A), and $C_{max}$ is calculated from $\langle C \rangle$ using Eq. 17. We emphasize that all three approaches account for the possibility that the measured cell is the product of cell division that occurred during the SIP incubation.

The assumption underlying Approaches A and B is that the measured cell assimilated C according to first-order kinetics (i.e., at a constant *k*) during the SIP incubation. Because $\langle C \rangle$ in Eq. 21 is the average C content of the measured cell across the cell cycle, *r* calculated by Approach A represents the *average* cell-specific rate across the cell cycle. Thus, inherent to Approach A is the assumption that the measured cell assimilated C according to first-order kinetics across the entire cell cycle and with the same value of *k* as determined from the SIP experiment (i.e., by Eq. 20). In contrast, *r* calculated by Approach B represents the *instantaneous* cell-specific C assimilation rate at the time of sampling. Note that the assumption of first-order kinetics implies that *r* would have been changing during the SIP experiment (because $r = k \cdot C$, where *k* is constant, and *C* is time-dependent; Figure 2E). Thus, *r* calculated by Approach B represents a rate at a specific time point during the cell cycle of the measured cell, namely at the end of the SIP experiment, when the C content reached the value of *C* used in Eq. 22.

The assumption underlying Approach C is that the measured cell assimilated C according to zero-order kinetics (i.e., at a constant *r*) during the SIP incubation and that the C content across the cell cycle varies strictly between the critical values of $C_{max}/2$ and $C_{max}$. Because this approach accounts for the possibility that the measured cell is the product of cell division that occurred during the SIP incubation, *r* calculated by Approach C represents the average cell-specific rate during the SIP incubation. Note that because Approach C makes use of both *C* and $\langle C \rangle$, it allows the reconstruction of the cell cycle stage of the measured cell. Thus, if the SIP experiment is conducted over a time interval that is short compared with the doubling time of cells in a population, approach C can potentially reveal how the cell-specific rates vary across the cell cycle.

## Application to Hypothetical Data

To illustrate the utility of the procedure described above, we calculated C assimilation rates based on hypothetical data from
## SIP experiment + nanoSIMS measurement

<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

$x_{S,eff}$: ¹³C atom fraction of the effective C source

## What is known about C content of the measured cell?

⟨C⟩: average C content across the full cell cycle§

C: C content at the end of SIP experiment

§ Equal to the average C content in a population with perfectly unsynchn-  
ronised cell cycles. It can also be calculated from the cellular C content  
at the time of cell division (Cmax) using the following equations:

* Zero-order kinetics of C assimilation (constant r):

$$
\langle C \rangle = C_{\text{max}} / (2 \times \ln(2)) \rightarrow C_{\text{max}} \approx 1.386 \times \langle C \rangle
$$

* First-order kinetics of C assimilation (constant k):

$$
\langle C \rangle = C_{\text{max}} \times \ln(2) \rightarrow C_{\text{max}} \approx 1.443 \times \langle C \rangle
$$

Calculate the C assimilation rate

k: Carbon-specific rate [h⁻¹]

r: Cell-specific rate [mol C h⁻¹ cell⁻¹]

Step 1

$$
x_S^E = \frac{x - x_i}{x_{S,eff} - x_i}
$$

Step 2

$$
k = -\frac{1}{t} \cdot \ln[1 - x_S^E]
$$

Step 3A

Only ⟨C⟩ can be constrained

$$
r = k \cdot \langle C \rangle
$$

r = average rate across cell cycle

Step 3B

Only C can be constrained

$$
r = k \cdot C
$$

r = instantaneous rate at the end of SIP experiment

Step 3C

Both ⟨C⟩ and C can be constrained

$$
r = Z^{-1}(t, C, C_{\text{max}}, x_S^E)
$$

r = average rate during SIP experiment

Underlying assumptions:

- All steps: Constant C density during the SIP experiment.

- Step 3A: *k* is constant across the cell cycle.

- Step 3C: Cell assimilated C according to zero-order kinetics, *r* was constant during the SIP experiment, cell divided into two identical cells when its C content reached *C*<sub>max</sub>

- Steps 2, 3A, 3B: Cell assimilated C according to first-order kinetics.

- Step 3B: *k* was constant during the SIP experiment.

**FIGURE 5** | Flowchart for quantifying C assimilation rates in single cells. First, cells are incubated with a ¹³C-labeled substrate and their ¹³C atom fractions, x, are measured by nanoSIMS. Second, the C content of the cells is constrained. Rate calculation then proceeds in three steps: (1) convert x to x<sub>S</sub><sup>E</sup> to account for the isotopic composition of the substrate, (2) calculate the carbon-specific C assimilation rate, k (this rate is independent of the C content of the cell), and (3) calculate the cell-specific C assimilation rate, r. Depending on how well the C content of the measured cell can be constrained, the third step proceeds either by using the value of k from step 2 (Approach A and B) or using the value of x<sub>S</sub><sup>E</sup> from step 1 and the inverse of the zig-zag function, Z<sup>-1</sup> (Approach C). The calculated r then represents the average rate across the cell cycle (Approach A), the instantaneous rate at the end of the SIP experiment (Approach B), or the average rate during the SIP experiment (Approach C). Steps 2 and 3 consider that the measured cell may be a product of a cell that divided during the incubation. The calculation approaches are implemented in a Matlab script available via GitHub (https://github.com/Ipolerecky/LARS).

cells incubated with a ¹³C-labeled C source ($x_{S,eff}$ = 0.1) for 2 hours. Measured cells had low, intermediate, and high levels of ¹³C labeling ($x_1$ = 0.02, $x_2$ = 0.033, $x_3$ = 0.09). Using Steps 1 and 2, these values yield carbon-specific assimilation rates of $k_1$ = 0.053 $h^{-1}$, $k_2$ = 0.14 $h^{-1}$, and $k_3$ = 1.1 $h^{-1}$, which correspond to doubling times of about $\tau_1$ = 13 h, $\tau_2$ = 4.9 h, and $\tau_3$ = 0.63 h, respectively.

To calculate the cell-specific rate, *r*, we further assumed three values for the C content of the measured cell (C = 7.6, 10.4 and 13.2 fmol), each quantified with the same analytical precision of ΔC = 0.1 fmol, and an average C content of ⟨C⟩ = 10 fmol (which corresponds to the critical C content of Cmax = 13.9 fmol; Eq. 17). **Figure 6** shows that *r* values predicted using Approach A do not depend on the C content of the measured cell (black bars), whereas they increase proportionally to the C content when

calculated using Approach B (gray bars). This pattern is expected from Eqs. 21 and 22. In contrast, *r* calculated by Approach C (open bars) does not necessarily vary monotonously with the C content of the measured cell. This occurs if the number of divisions during the SIP incubation differs between cells with the same ¹³C labeling but with different C contents at the end of the incubation. For example, the cell with C = 7.6 fmol divided once, whereas the cell with C = 10.4 fmol did not divide during the SIP incubation (Figure 6B, open bars). Thus, the first cell needs to assimilate C more rapidly than the second cell so that both reach the same ¹³C labeling at the end of the incubation. Differences between the *r* values calculated by Approach A and C tend to diminish with increasing ¹³C labeling of the cells (compare black and open bars in Figures 6A–C). This pattern emerges because the difference between the first-order and zero-order models
of C assimilation decreases with the number of cell divisions during the incubation and thus with the increasing extent of ¹³C labeling (Figure 3).

## DISCUSSION

This work uses mathematical modeling to identify and explicitly formulate fundamental assumptions underpinning calculations of substrate assimilation rates into cells based on SIP-nanoSIMS data. The 3-step approach to calculate substrate assimilation rates presented in this research evaluates the impact of the model of C assimilation by a cell (zero-order or first-order kinetics), cell division, and the carbon content of the cell measured by nanoSIMS on the calculated rate. Below, we summarize the key assumptions underlying this 3-step approach and compare it to previously published approaches for determining rates of substrate assimilation. Then, we discuss how the rates evaluated in single cells can be upscaled toward a cell population. Finally, motivated by results from a recent nanoSIMS study on diazotrophic cyanobacteria (Polerecky et al., 2021), we identify key sources for the variation in calculated assimilation rates typically observed among cells of the same species within a population that have not been previously discussed.

## Summary of Assumptions Used to Calculate Assimilation Rates

Important assumptions underlying the calculations described in this study are: (1) Catabolism or turnover of cell material, such as protein and carbohydrates, does not occur even though their half-lives may be shorter than the incubation duration. Thus, the cell-specific C assimilation rate directly reflects the increase in the cell biomass over time (Eq. 1). (2) The cell's activity is constant during the incubation. That is, the C assimilation rates *r* and *k* in the zero-order and first-order kinetic models, respectively, are time-independent. (3) Carbon density of the cells is constant or unaffected by accumulation of storage inclusions. (4) A cell divides into two identical cells when its C content reaches a certain critical value *C*<sub>max</sub>. (5) The effect of kinetic isotope fractionation is minor and negligible. (6) Isotopic composition of the substrate, *x*<sub>S,eff</sub>, does not change during the incubation. Analysis of a SIP incubation with time-dependent *x*<sub>S,eff</sub> is discussed in more detail in Supplementary Methods, Section "Modeling Cellular Assimilation of Substrates With Time-Dependent Isotopic Labeling."

## Comparison With Previous Studies

Commonly used approaches for calculating substrate assimilation rates (Foster et al., 2013; Krupke et al., 2015; Svedén et al., 2015; Klawonn et al., 2016; Stryhanyuk et al., 2018; Calabrese et al., 2019; Mills et al., 2020) assume a zero-order kinetic model of substrate assimilation and do not consider that the measured cell could be a product of cell division during the SIP incubation. Additionally, they typically use an average value for the C content or the C density combined with the biovolume of the measured cell. Here we evaluate the impact of these assumptions on calculated rates.

The most common published approaches calculate *r* according to formulas equivalent to (Krupke et al., 2015; Svedén et al., 2015; Klawonn et al., 2016; Mills et al., 2020; Ploug, 2021)

$$r = \frac{1}{t} \cdot x_S^E \cdot C_f, \qquad (24)$$

or (Stryhanyuk et al., 2018; Calabrese et al., 2019)

$$r = \frac{1}{t} \cdot \frac{x_S^E}{1-x_S^E} \cdot C_i. \qquad (25)$$

In these formulas, *C*<sub>i</sub> and *C*<sub>f</sub> denote the C content of the measured cell at the beginning and end of the SIP incubation, respectively. Equations 24 and 25 are equivalent and follow directly from the mass balance of a non-dividing cell. Equation 24 is obtained by defining the cell-specific C assimilation rate as *r* = *C*<sub>a</sub>/*t*, where *C*<sub>a</sub> is the amount of carbon assimilated by the cell during time *t*, and by considering that *x*<sub>S</sub><sup>E</sup> = *C*<sub>a</sub>/*C*<sub>f</sub> for a cell that did not divide during time *t* (see Supplementary Methods, Section "Relationships between x<sub>S</sub><sup>E</sup>, X<sub>net</sub>, F<sub>Xnet</sub> and K<sub>A</sub>"). Similarly, Eq. 25 derives from the same assumptions and by additionally considering that, for a cell that did not divide, *C*<sub>i</sub> and *C*<sub>f</sub> are related according to *C*<sub>f</sub> = *C*<sub>i</sub> + *C*<sub>a</sub>. In summary, these equations yield the average cell-specific C assimilation rate during the SIP incubation. However, the calculated *r* is accurate only if the cell did not divide during the incubation and if *C*<sub>f</sub> and *C*<sub>i</sub> are stringently constrained.

Constraining *C*<sub>i</sub> or *C*<sub>f</sub> of a cell is not trivial even without considering cell division. Facing the practical difficulties of their direct measurement, some practitioners use ⟨C⟩ in place of *C*<sub>f</sub> (e.g., Svedén et al., 2015; Olofsson et al., 2019a) or approximate *C*<sub>i</sub> using *C*<sub>f</sub> (Stryhanyuk et al., 2018; Calabrese et al., 2019). These substitutions, however, invalidate the approach based on Eq. 24 or 25 because the mass balance underpinning these formulas is no longer achieved.

Many studies estimate the biovolume of the measured cells from nanoSIMS images and apply literature values for the C density, or an empirical relationship between the C density and biovolume of phytoplankton cells (e.g., Verity et al., 1992; Stryhanyuk et al., 2018; Khachikyan et al., 2019), to approximate *C*<sub>f</sub> (e.g., Foster et al., 2013; Krupke et al., 2015; Schoffelen et al., 2018; Mills et al., 2020; Trembath-Reichert et al., 2021). Although this approach may yield a well-constrained value of *C*<sub>f</sub>, one must be careful when using it in Eq. 24 to calculate *r*. Specifically, if the cell divided during the SIP incubation, the correct value of *C*<sub>f</sub> should include the C content of the *measured* cell as well as its *sister* cell (or sister cells, if the cell divided more than once). Only then will the mass-balance requirement be fulfilled, and Eq. 24 will yield the correct value of *r* for the measured cell. This could, in principle, be achieved by including a correction factor when calculating the *cell*-specific rate using Eq. 24. Several recent studies (Schoffelen et al., 2018, 2019; Olofsson et al., 2019a,b; Ploug, 2021) have employed a factor 2 in the calculation of *carbon*-specific growth rates. However, as this latter quantity is solely determined by the cell's isotopic composition (i.e., it does not describe absolute amounts of carbon assimilated; see Eq. 20),
<figure><img src="image_11.png" /></figure>

**FIGURE 6** | Cell-specific rates of C assimilation in cells determined using hypothetical data. The rates were calculated using the procedure described in the text (Steps 1-3) and using approaches A, B and C for Step 3. Each panel (A–C) corresponds to a different value of the ¹³C atom fraction of the cell, *x*, as shown above the graph. The rates represent the average rate across the cell cycle (Approach A), the instantaneous rate at the end of the SIP incubation (Approach B), and the average rate during the SIP incubation (Approach C). Values were calculated for (C) = 10 fmol C cell⁻¹, xS,eff = 0.1, *t* = 2 h, and three values of the C content of the measured cell (x-axis), each determined with the precision of 0.1 fmol. Numbers in brackets indicate how many times the measured cell divided during the incubation, as determined by the algorithm implementing Approach C. Note differences in the y-scale.

this operation does not correct for the carbon content of a sister cell as described above.

To illustrate quantitatively the impact of cell division on calculated *r*, we compared the results obtained using Eq. 24 and our Approach C (Eq. 23) for a range of ¹³C labeling and C content of the measured cell (Figure 7). We did not apply a correction factor to calculate *C*<sub>f</sub> but assumed that *C*<sub>f</sub> in Eq. 24 is the same as C in Eq. 23. As expected, both approaches yield the same result if the cell did not divide during the SIP incubation (see values in Figure 7 for which the number of cell divisions is equal to 0), because the function in Eq. 11, which is the basis for Eq. 24, is equal to the zig-zag function Z (see Figure 2A, black curve before *t*/*τ* < 1 and red curve before *t*/*τ* < 0.5). However, if the measured cell did divide during the incubation, Eq. 24 yields a lower *r* value compared to that calculated by Approach C (Figure 7). The difference between the two approaches becomes more pronounced as *x*<sub>S</sub><sup>E</sup> of the measured cell increases above a threshold value that depends on the cell cycle stage (Figure 7B). For example, the threshold value of *x*<sub>S</sub><sup>E</sup> is 0.1 if the cell cycle stage of the measured cell is 10%, but it increases to 0.33 and 0.47 for the cell cycle stage of 50% and 90%, respectively (inset in Figure 7B).

The question remains, however, how to determine whether a cell measured by nanoSIMS has divided during the SIP incubation. Tracing the history of a cell analyzed by nanoSIMS is practically impossible because nanoSIMS is a terminal measurement. The history can be reconstructed if additional information about the measured cell, namely ⟨C⟩, can be constrained. Constraining ⟨C⟩ allows evaluation of the minimum and maximum C content of the measured cell (C<sub>max</sub>/2 and C<sub>max</sub>, respectively; Eq. 17). Combining this information with the carbon content, C, enables estimation of the cell cycle stage of the measured cell. Combined with the value of x<sub>S</sub><sup>E</sup>, the history of the measured cell during the SIP experiment can be reconstructed to estimate the value of *r* that accounts for cell division. This sequence of steps is implemented by our Approach C based on the inverse of the zig-zag function (Eq. 23).

We conclude that the 3-step procedure proposed in this study (Figure 5) is applicable to many situations encountered by

users of SIP-nanoSIMS to estimate substrate assimilation rates in individual cells. When the C content of a measured cell, *C*, cannot be constrained, the carbon-specific assimilation rate, *k* (Eq. 20), and the corresponding doubling time, τ = (ln2)/*k*, is the maximal information that can be estimated from the measured isotopic composition of the cell. Constraining *C* allows estimation of the instantaneous cell-specific assimilation rate at the end of the SIP experiment (Approach B, Eq. 22), while constraining ⟨C⟩ allows extrapolation towards the average cell-specific assimilation rate across the cell cycle (Approach A, Eq. 21). Finally, constraining both *C* and ⟨C⟩ allows determination of the average cell-specific assimilation rate over the SIP experiment (Approach C, Eq. 23). Note that accounting for cell division, as done in Approach C, critically depends on the assumption that the measured cell divided during the incubation when its *C* content reached *C*<sub>max</sub> estimated from ⟨C⟩. If this assumption is not valid (e.g., because the *C* assimilated by the cell was allocated into storage inclusions, or cell division was delayed or not binary), or cannot be verified, the instantaneous rate at the end of the SIP experiment (Eq. 22) is the maximal information about the cell-specific assimilation rate that can be estimated by combining the ¹³C labeling and *C* content of the measured cell.

## Extrapolation Toward a Cell Population

A frequent aim of SIP-nanoSIMS measurements is to estimate a bulk C assimilation rate for a target cell population and thus evaluate its potential impact on C fluxes in the environment. This upscaling can be done by first averaging the cell-specific C assimilation rates determined for individual target cells to estimate the population average, ⟨r⟩, and then multiplying ⟨r⟩ by the total target cell abundance in the environment. When calculating the assimilation rates for individual cells, Approach A or B is recommended if cell cycles in the target population are perfectly unsynchronized, whereas Approach C is recommended for partially synchronized cells. Because these approaches account for cell division, ⟨r⟩ estimated in this way will reflect the population average across a broad range of incubation times (Supplementary Figure 5), which may be important when using SIP-nanoSIMS to study environmental
<figure><img src="image_12.png" /><figcaption>**FIGURE 7** | Comparison of previously published and the current approaches for calculating C assimilation rates in cells. The key difference is that our Approach C accounts for cell division, whereas cell division is not considered in the approach using Eq. 24. (A) Carbon-specific rates predicted from the ¹³C labeling of the measured cell, expressed as $x_S^E$. Calculations assumed $\langle C \rangle = 10 \text{ fmol C cell}^{-1}$, $t = 2 \text{ h}$, and three values of the cell's C content, C. The values of C = 7.6, 10.4, and 13.2 fmol correspond to the cell cycle stage of s = 0.1, 0.5, and 0.9, respectively. (B) Ratio between the r value predicted by Approach C and the approach based on Eq. 24, calculated for several values describing the cell cycle stage of the measured cell (indicated by numbers next to the lines). Approach based on Eq. 24 yields lower rates compared to Approach C if the measured cell divided during the incubation [see number of divisions in the top graph in panel (A)]. This occurs when $x_S^E$ exceeds a threshold value that increases with the increasing cell cycle stage of the measured cell (see inset in panel (B)). Note that the ratio only depends on $x_S^E$ and s and not on the incubation time.</figcaption></figure>

samples where the growth rates of cells are *a-priori* unknown. Averaging of cell-specific rates calculated using Eq. 24 is not recommended if the expectation is that a significant number of target cells divided during the incubation, because this approach would underestimate $\langle r \rangle$ (Supplementary Methods, Section “Simulating SIP Incubation,” Supplementary Figure 5). Note that the precision and accuracy of the estimated bulk assimilation rate will still depend on the number of target cells measured by nanoSIMS, i.e., on the quality of sampling of the target population.

## Interpretation of Calculated Rates and Their Variability Among Cells

Typically, ¹³C atom fractions obtained by nanoSIMS measurements vary among cells, raising the question, to what extent this variability is caused by differences in the intrinsic metabolic activities of cells probed by the SIP experiment (i.e., *r*) versus differences arising from cells at different stages of their cell cycle during the labeling interval. Additionally, we note that assimilation rates obtained in nanoSIMS studies can reflect variation not only in the targeted process (e.g., C fixation) but also in other metabolic processes and cellular characteristics (e.g., simultaneous assimilation of other external or internal C sources). In this study we identified that in addition to *r*, the sources of variability include *t*, $x_{S,eff}$, and $C_i$. Our 3-step procedure (Figure 5) accounts for some, but not all, of these sources of variability. Below, we discuss what this uncertainty implies for interpreting the calculated rates of C assimilation (*k* or *r*) and especially their variability among measured cells. We focus on C, but similar arguments can be made for other elements.

### Influence of an Alternative C Source

The effective isotopic composition of the C source, $x_{S,eff}$, will almost certainly be different from the isotopic composition of the target C source, $x_S$, if cells assimilate carbon from an alternative source in addition to the target source. The parameter $f_{tar}$ accounts for the influence of alternative C sources (see Materials and Methods and Table 1). Calculated *k* values will underestimate total C assimilation if an alternative C source influenced $x_{S,eff}$ (Table 2). Thus, cell-to-cell differences in *k* may reflect intercellular variation in use of the target source rather than total C assimilation. Cell-specific rates, *r*, will be similarly affected as they are related to *k*.

### Influence of C Storage Content

Because *k* is defined as $k = r/C$, variability in calculated *k* among cells can be caused not only by variation in the instantaneous C assimilation rate (*r*) but also by variation in the instantaneous C content (*C*). Such variations in *C* can be expected, for example, in cells that encounter starvation, i.e., they have depleted a required growth resource from the environment, or the environment has

<table><caption>**TABLE 2** | Simulated effect of C assimilation from two sources.</caption><thead><tr><th rowspan="2">cell</th><th rowspan="2">f<sub>tar</sub></th><th rowspan="2">x</th><th colspan="2">x<sub>S,eff</sub></th><th colspan="2">k (h<sup>-1</sup>)</th></tr><tr><th>Correct</th><th>Assumed in calculation</th><th>Calculated</th><th>Total</th></tr></thead><tbody><tr><td>1</td><td>0.5</td><td>0.0197</td><td>0.505</td><td>0.99</td><td>0.005</td><td>0.01</td></tr><tr><td>2</td><td>0.8</td><td>0.0255</td><td>0.802</td><td>0.99</td><td>0.008</td><td>0.01</td></tr><tr><td>1</td><td>0.5</td><td>0.0197</td><td>0.505</td><td>0.505</td><td><strong>0.01</strong></td><td><strong>0.01</strong></td></tr><tr><td>2</td><td>0.8</td><td>0.0255</td><td>0.802</td><td>0.802</td><td><strong>0.01</strong></td><td><strong>0.01</strong></td></tr></tbody></table>

Shown are results for two cells that assimilate C with equal total carbon-specific rates of 0.01 h⁻¹ over 2 h. The cells differ only in the fraction of C assimilated from the target C source, $f_{tar}$. If C is incorrectly assumed to be assimilated from only the target source, calculated *k* will underestimate the total C assimilation rate (as discussed in Section Influence of an Alternative C Source). The calculated and total values will only be equal if C assimilation from the alternative C source is accounted for by using the correct value for $x_{S,eff}$ (lines 3-4, values shown in bold). Calculated for $x_{S,tar} = 0.99$, $x_{S,alt} = 0.011$, and $x_i = 0.011$.
changed, such that assimilation of that resource is no longer
possible. These conditions cause cells to cease exponential growth
and commonly lead to the synthesis of storage compounds (e.g.,
polysaccharides). The variability in calculated *k* (Eq. 20) among
starving cells may therefore reflect variation in storage content
rather than *r*.

Similarly, cell-to-cell variability in calculated *r* needs to be
interpreted with caution as it may reflect variation in storage
content, rather than intrinsic cell-specific C assimilation rates,
because neither Approach A, B nor C accounts for potential
variation in the cell's C content caused by storage (see the
underlying assumptions in **Figure 5**). Specifically, *r* is calculated
either directly from *k* using ⟨C⟩ (Approach A, Eq. 21), or the
calculation relies on the assumption that the cell C density, which
may be used to estimate the C content of the measured cell based
on its biovolume (**Supplementary Methods**, Section "Estimating
Cellular C Content"), is constant during the SIP experiment
and equal among cells (Approach B and C, Eq. 22-23). These
assumptions are unlikely to be valid in cells that metabolize
storage inclusions. By accounting for the relationships between
storage content and apparent assimilation rates, however, we can
deduce information about storage metabolism from nanoSIMS
data, as outlined in the following section.

Interpretation of Simultaneous <sup>13</sup>C and <sup>15</sup>N Labeling Experiments

For a single labeled element (e.g., ¹³C or ¹⁵N), the importance of storage metabolism versus assimilation of new material is difficult to resolve without additional measurements at subcellular level (e.g., TEM images of storage inclusions). On the other hand, dual labeling (i.e., simultaneous incubations with two labeled substrates such as ¹³C and ¹⁵N) can greatly increase informative value of the SIP experiment. By systematic consideration of the effects of storage inclusion metabolism on the element-specific rates of C and N assimilation (*k*<sub>C</sub> and *k*<sub>N</sub>), whole-cell nanoSIMS measurements can be used to gain information not only about the acquisition of new carbon and nitrogen, but also about intracellular fluxes of these elements (Polerecky *et al.*, 2021). Specifically, metabolic pathways leading to storage inclusion biosynthesis and catabolism can be revealed by leveraging information about the known elemental compositions of specific storage inclusions (e.g., polysaccharides or cyanophycin) and the isotopic compositions of target C and N sources in carefully designed SIP experiments.

We simulated assimilation of ¹³C- and ¹⁵N-labeled substrates by cells under various scenarios, including preferential incorporation of the assimilated C and N into different cell compartments (e.g., cell matrix or storage inclusions), or C and N incorporation using internally recycled C and N in addition to C and N originating from the target sources provided externally (Supplementary Methods, Section “Modeling Simultaneous Assimilation of C and N”). Subsequently, we calculated the C- and N-specific assimilation rates, *k*<sub>C</sub> and *k*<sub>N</sub>, in these cells using the approach described above (Eq. 20).

If a cell assimilated the ¹³C- and ¹⁵N-labeled source substrates
in a balanced way, i.e., the amounts of assimilated C and N did
not change the C:N ratio of the cell, the ratio of the calculated
element-specific rates of C and N assimilation, *k*<sub>C</sub>*k*<sub>N</sub>, is equal

to 1. Deviation of *k*<sub>C</sub>*k*<sub>N</sub> from 1 can arise for multiple reasons. First, the effective isotopic composition of the C or N source, x(¹³C)<sub>S,eff</sub> or x(¹⁵N)<sub>S,eff</sub>, differed from the isotopic composition of the target C and N source (Figure 8, stars). For example, by providing autotrophs with ¹³C-DIC (their only carbon source) and ¹⁵N-ammonium, a value of *k*<sub>C</sub>*k*<sub>N</sub> > 1 could indicate uptake of organic N compounds in addition to ammonium. Thus, deviation of *k*<sub>C</sub>*k*<sub>N</sub> from 1 can reveal uptake of additional substrates. A similar approach was used by Schoffelen et al. (2018) to reveal uptake of organic P in addition to inorganic P. Second, C and N assimilation did not match the C:N ratio of the whole cell. For example, ¹³C or ¹⁵N were incorporated into inclusions with C:N ratios that differed from that of the whole cell (Figure 8, circle, square, triangle). This can happen in cells that temporally decouple assimilation of C and N, such as in some diazotrophs, which assimilate C in the light and N mostly in the dark (e.g., Polerecky et al., 2021). Third, C or N from existing cell material

<figure><img src="image_13.png" /></figure>

**FIGURE 8** | Simulations of cell-specific assimilation of C and N. Element-specific rates of C and N assimilation by whole cells, *k*<sub>C</sub> and *k*<sub>N</sub>, calculated from hypothetical data. The first simulation (stars) shows the results for assimilated C and N that originated from two sources: labeled target sources [x(¹³C)<sub>S,tar</sub> = 1, x(¹⁵N)<sub>S,tar</sub> = 1] and unlabeled alternative sources [x(¹³C)<sub>S,alt</sub> = 0.011, x(¹⁵N)<sub>S,alt</sub> = 0.0037]. Resultant *k*<sub>C</sub> and *k*<sub>N</sub> values assumed assimilation of only the target sources [x(¹³C)<sub>S,eff</sub> = 1, x(¹⁵N)<sub>S,eff</sub> = 1]. Individual data-points correspond to cells that assimilated different fractions of C and N from the target sources, with corresponding *f*<sub>C,tar</sub> and *f*<sub>N,tar</sub> values given in parentheses (as *f*<sub>C,tar</sub>*f*<sub>N,tar</sub>). Deviation of *k*<sub>C</sub>*k*<sub>N</sub> from 1 reveals assimilation of C or N from alternative sources. The second simulation shows the impacts of preferential incorporation of C and N into cyanophycin granules (cy, circle), polysaccharide inclusions (ps, triangle), and cell matrix (m, square). Cellular compartments with C:N that differ from the whole cell can cause *k*<sub>C</sub>*k*<sub>N</sub> to deviate from 1. The third simulation (diamonds) shows results for cyanophycin synthesis using unlabeled C and N from inside the cell (e.g., protein and polysaccharide catabolism) in addition to C and N from the external target sources. The parentheses give *f*<sub>C,tar</sub>*f*<sub>N,tar</sub>. Use of unlabeled internal C and N pools causes the cell to assimilate C and N from sources with x(¹³C)<sub>S,eff</sub> < 1 and x(¹⁵N)<sub>S,eff</sub> < 1, whereas calculation of *k*<sub>C</sub> and *k*<sub>N</sub> assumed x(¹³C)<sub>S,eff</sub> = 1 and x(¹⁵N)<sub>S,eff</sub> = 1. Details of the simulations are given in **Supplementary Methods**, Section “Modeling Simultaneous Assimilation of C and N.”
(e.g., proteins, storage compounds) was mobilized into a new form. For example, protein degradation released N that was then used for cyanophycin synthesis in addition to incorporation of ¹⁵N from the target source (Figure 8, diamonds; see also Polerecky et al., 2021). Thus, deviation of $k_C/k_N$ from 1 can be used to identify the presence of intra-cellular C and N stores and evaluate their synthesis or mobilization without the need for direct measurement at the sub-cellular level. Notably, the assessment of the deviation of $k_C/k_N$ from 1 is solely based on the nanoSIMS measurement and does not require knowledge about the C and N content of the measured cell.

## CONCLUSION

This study revisited calculations of biologically relevant parameters from SIP-nanoSIMS data. We suggest a step-by-step procedure to evaluate elemental assimilation rates into whole cells, discuss key assumptions underlying the procedure, and describe factors that can introduce variability into the calculated rates.

Determining single-cell assimilation rates is challenging, as it requires proper evaluation of multiple factors that are difficult to constrain experimentally or analytically, necessitating assumptions in the calculation process. Our analysis illustrates and quantifies the consequences of some of the key assumptions, particularly focusing on assimilation kinetics (zero vs. first order), cell division, and the elemental content of the measured cell. We give guidance for interpreting the rates calculated by approaches that differ depending on knowledge of the elemental content of the measured cell and develop an approach to account for cell division in the calculation process.

We call for caution when interpreting inter-cellular variability of assimilation rates calculated from SIP-nanoSIMS data if (1) cells can assimilate elements such as C from alternative sources present in the environment in addition to the target source added externally during the SIP experiment, (2) storage inclusions are synthesized, catabolized or merely present inside the studied cells during the SIP experiment, or (3) cell divisions during the SIP experiment may be expected. Our modeling analysis further shows that for dual labeling SIP experiments (e.g., using ¹³C and ¹⁵N labeled substrates), deviation of $k_C/k_N$ from 1 can be used to identify the presence of intra-cellular C and N stores and evaluate their synthesis or mobilization without the need for direct measurement at the sub-cellular level.

## REFERENCES

* Ackermann, M. (2015). A functional perspective on phenotypic heterogeneity in microorganisms. *Nat. Rev. Microbiol.* 13, 497–508. doi: 10.1038/nrmicro3491

* Arandia-Gorostidi, N., Weber, P. K., Alonso-Sáez, L., Morán, X. A. G., and Mayali, X. (2017). Elevated temperature increases carbon and nitrogen fluxes between phytoplankton and heterotrophic bacteria through physical attachment. *ISME J.* 11, 641–650. doi: 10.1038/ismej.2016.156

* Berthelot, H., Duhamel, S., L'Helguen, S., Maguer, J.-F., Wang, S., Cetinić, I., et al. (2019). NanoSIMS single cell analyses reveal the contrasting nitrogen sources for small phytoplankton. *ISME J.* 13, 651–662. doi: 10.1038/s41396-018-0285-8

## DATA AVAILABILITY STATEMENT

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author/s.

## AUTHOR CONTRIBUTIONS

LP conceived the study and performed the mathematical analysis and modeling. LP, KH, and ME wrote the manuscript with significant contributions from all other co-authors. All authors contributed with ideas that shaped the final version of the manuscript.

## FUNDING

ME was supported by the Czech Science Foundation (GA CR, Grant Number 20-02827Y). TM was supported by the Czech Science Foundation (GA CR, Grant No 20-17627S). TZ was supported by the Ministry of Education, Youth and Sports of the Czech Republic (OP RDE Grant No. CZ.02.1.01/0.0/0.0/16-026/0008413 'Strategic Partnership for Environmental Technologies and Energy Production') and by the Czech Science Foundation (GA CR, Grant No. 18-24397S). DC was supported by Mobility project: CZ.02.2.69/0.0/0.0/16_027/0007990, International mobility of researchers of the Institute of Microbiology of the CAS.

## ACKNOWLEDGMENTS

This work was inspired by discussions and observations that emerged from the Group for Aquatic Photosynthesis (GAP) workshop held in Třeboň, Czech Republic, in August 2017. We are particularly grateful to Ondřej Prášil for organizing the workshop. We also thank Celine Loussert-Fonta, Marcel Kuypers, Peter K. Weber, and Martin Von Bergen whose constructive comments helped improve the manuscript.

## SUPPLEMENTARY MATERIAL

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb.2021.621634/full#supplementary-material

* Bonnet, S., Berthelot, H., Turk-Kubo, K., Cornet-Barthaux, V., Fawcett, S., Berman-Frank, I., et al. (2016). Diazotroph derived nitrogen supports diatom growth in the South West Pacific: a quantitative study using nanoSIMS. *Limnol. Oceanogr.* 61, 1549–1562. doi: 10.1002/lno.10300

* Boschker, H. T. S., and Middelburg, J. J. (2002). Stable isotopes and biomarkers in microbial ecology. *FEMS Microbiol. Ecol.* 40, 85–95. doi: 10.1111/j.1574-6941.2002.tb00940.x

* Calabrese, F., Voloshynovska, I., Musat, F., Thullner, M., Schlömann, M., Richnow, H. H., et al. (2019). Quantitation and comparison of phenotypic heterogeneity among single cells of monoclonal microbial populations. *Front. Microbiol.* 10:2814. doi: 10.3389/fmicb.2019.02814
Collins, J. F., and Richmond, M. H. (1962). Rate of growth of *Bacillus cereus* between divisions. *J. Gen. Microbiol.* 28, 15-33.

Coplen, T. B. (2011). Guidelines and recommended terms for expression of stable-isotope-ratio and gas-ratio measurement results. *Rapid Commun. Mass Spectrom.* 25, 2538-2560. doi: 10.1002/rcm.5129

Dekas, A. E., Parada, A. E., Mayali, X., Fuhrman, J. A., Wollard, J., Weber, P. K., et al. (2019). Characterizing chemoautotrophy and heterotrophy in marine Archaea and bacteria with single-cell multi-isotope NanoSIP. *Front. Microbiol.* 10:2682. doi: 10.3389/fmicb.2019.02682

Dumont, M. G., and Murrell, J. C. (2005). Stable isotope probing—linking microbial identity to function. *Nat. Rev. Microbiol.* 3, 499–504. doi: 10.1038/nrmicro1162

Eichner, M. J., Klawonn, I., Wilson, S. T., Littmann, S., Whitehouse, M. J., Church, M. J., et al. (2017). Chemical microenvironments and single-cell carbon and nitrogen uptake in field-collected colonies of *Trichodesmium* under different pCO2. *ISME J.* 11, 1305–1317. doi: 10.1038/ismej.2017.15

Finzi-Hart, J. A., Pett-Ridge, J., Weber, P. K., Popa, R., Fallon, S. J., Gunderson, T., et al. (2009). Fixation and fate of C and N in the cyanobacterium *Trichodesmium* using nanometer-scale secondary ion mass spectrometry. *Proc. Natl. Acad. Sci. U.S.A.* 106, 6345–6350. doi: 10.1073/pnas.0810547106

Foster, R. A., Kuypers, M. M. M., Vagner, T., Paerl, R. W., Musat, N., and Zehr, J. P. (2011). Nitrogen fixation and transfer in open ocean diatom–cyanobacterial symbioses. *ISME J.* 5, 1484–1493. doi: 10.1038/ismej.2011.26

Foster, R. A., Sztejrenszus, S., and Kuypers, M. M. M. (2013). Measuring carbon and N2 fixation in field populations of colonial and free-living unicellular cyanobacteria using nanometer-scale secondary ion mass spectrometry. *J. Phycol.* 49, 502–516. doi: 10.1111/jpy.12057

Geerlings, N. M. J., Geelhoed, J. S., Vasquez-Cardenas, D., Kienhuis, M. V. M., Hidalgo-Martinez, S., Boschker, H. T. S., et al. (2021). Cell cycle, filament growth and synchronized cell division in multicellular cable bacteria. *Front. Microbiol.* 12:620807. doi: 10.3389/fmicb.2021.620807

Geerlings, N. M. J., Karman, C., Trashin, S., As, K. S., Kienhuis, M. V. M., Hidalgo-Martinez, S., et al. (2020). Division of labor and growth during electrical cooperation in multicellular cable bacteria. *Proc. Natl. Acad. Sci. U.S.A.* 117:5478. doi: 10.1073/pnas.1916244117

Harding, K., Turk-Kubo, K. A., Sipler, R. E., Mills, M. M., Bronk, D. A., and Zehr, J. P. (2018). Symbiotic unicellular cyanobacteria fix nitrogen in the Arctic Ocean. *Proc. Natl. Acad. Sci. U.S.A.* 115:13371. doi: 10.1073/pnas.181365 8115

Hoppe, P., Cohen, S., and Meibom, A. (2013). NanoSIMS: technical aspects and applications in cosmochemistry and biological geochemistry. *Geoostand. Geoanal. Res.* 37, 111–154. doi: 10.1111/j.1751-908X.2013.00239.x

Hungate, B. A., Mau, R. L., Schwartz, E., Caporaso, J. G., Dijkstra, P., van Gestel, N., et al. (2015). Quantitative microbial ecology through stable isotope probing. *Appl. Environ. Microbiol.* 81, 7570–7581. doi: 10.1128/AEM.02280-15

Jehmlich, N., Vogt, C., Lünsmann, V., Richnow, H. H., and von Bergen, M. (2016). Protein-SIP in environmental studies. *Curr. Opin. Biotechnol.* 41, 26–33. doi: 10.1016/j.copbio.2016.04.010

Khachikyan, A., Milucka, J., Littmann, S., Ahmerkamp, S., Meador, T., Könneke, M., et al. (2019). Direct cell mass measurements expand the role of small microorganisms in nature. *Appl. Environ. Microbiol.* 85:e00493-19. doi: 10. 1128/AEM.00493-19

Klawonn, I., Nahar, N., Walve, J., Andersson, B., Olofsson, M., Svedén, J. B., et al. (2016). Cell-specific nitrogen- and carbon-fixation of cyanobacteria in a temperate marine system (Baltic Sea). *Environ. Microbiol.* 18, 4596–4609. doi: 10.1111/1462-2920.13557

Koch, A. L. (1966). Distribution of cell size in growing cultures of bacteria and the applicability of the collins-richmond principle. *J. Gen. Microbiol.* 45, 409–417. doi: 10.1099/00221287-45-3-409

Koch, A. L., and Schaechter, M. (1962). A model for statistics of the cell division process. *J. Gen. Microbiol.* 29, 435–454. doi: 10.1099/00221287-29-3-435

Krupke, A., Mohr, W., LaRoche, J., Fuchs, B. M., Amann, R. I., and Kuypers, M. M. M. (2015). The effect of nutrients on carbon and nitrogen fixation by the UCYN-A–haptophyte symbiosis. *ISME J.* 9, 1635–1647. doi: 10.1038/ismej.2014.253

Loussert-Fonta, C., Toullec, G., Paraecattil, A. A., Jeangros, Q., Krueger, T., Escrig, S., et al. (2020). Correlation of fluorescence microscopy, electron microscopy,

and NanoSIMS stable isotope imaging on a single tissue section. *Commun. Biol.* 3:362. doi: 10.1038/s42003-020-1095-x

Martínez-Pérez, C., Mohr, W., Löscher, C. R., Dekaezemacker, J., Littmann, S., Yilmaz, P., et al. (2016). The small unicellular diazotrophic symbiont, UCYN-A, is a key player in the marine nitrogen cycle. *Nat. Microbiol.* 1:16163. doi: 10.1038/nmicrobiol.2016.163

Matantseva, O., Skarlato, S., Vogts, A., Pozdnyakov, I., Liskow, I., Schubert, H., et al. (2016). Superposition of individual activities: urea-mediated suppression of nitrate uptake in the dinoflagellate prorocentrum minimum revealed at the population and single-cell levels. *Front. Microbiol.* 7:1310. doi: 10.3389/fmicb. 2016.01310

Mayali, X. (2020). NanoSIMS: microscale quantification of biogeochemical activity with large-scale impacts. *Ann. Rev. Mar. Sci.* 12, 449–467. doi: 10.1146/annurev- marine-010419-010714

Mills, M. M., Turk-Kubo, K. A., van Dijken, G. L., Henke, B. A., Harding, K., Wilson, S. T., et al. (2020). Unusual marine cyanobacteria/haptophyte symbiosis relies on N2 fixation even in N-rich environments. *ISME J.* 14, 2395–2406. doi: 10.1038/s41396-020-0691-6

Milucka, J., Ferdelman, T. G., Polerecky, L., Franzke, D., Wegener, G., Schmid, M., et al. (2012). Zero-valent sulphur is a key intermediate in marine methane oxidation. *Nature* 491, 541–546. doi: 10.1038/nature11656

Musat, N., Foster, R., Vagner, T., Adam, B., and Kuypers, M. M. M. (2012). Detecting metabolic activities in single cells, with emphasis on nanoSIMS. *FEMS Microbiol. Rev.* 36, 486–511. doi: 10.1111/j.1574-6976.2011.00300.x

Musat, N., Halm, H., Winterholler, B., Hoppe, P., Peduzzi, S., Hillion, F., et al. (2008). A single-cell view on the ecophysiology of anaerobic phototrophic bacteria. *Proc. Natl. Acad. Sci. U.S.A.* 105, 17861–17866. doi: 10.1073/pnas. 0809329105

Neufeld, J. D., Dumont, M. G., Vohra, J., and Murrell, J. C. (2007). Methodological considerations for the use of stable isotope probing in microbial ecology. *Microb. Ecol.* 53, 435–442.

Nuñez, J., Renslow, R., Cliff, J. B., and Anderton, C. R. (2017). NanoSIMS for biological applications: current practices and analyses. *Biointerphases* 13:03B301. doi: 10.1116/1.4993628

Olofsson, M., Kourtchenko, O., Zetsche, E. M., Marchant, H. K., Whitehouse, M. J., Godhe, A., et al. (2019a). High single-cell diversity in carbon and nitrogen assimilations by a chain-forming diatom across a century. *Environ. Microbiol.* 21, 142–151. doi: 10.1111/1462-2920.14434

Olofsson, M., Robertson, E. K., Edler, L., Arneborg, L., Whitehouse, M. J., and Ploug, H. (2019b). Nitrate and ammonium fluxes to diatoms and dinoflagellates at a single cell level in mixed field communities in the sea. *Sci. Rep.* 9:1424. doi: 10.1038/s41598-018-38059-4

Pett-Ridge, J., and Weber, P. K. (2012). “NanoSIP: NanoSIMS applications for microbial biology,” in *Microbial Systems Biology: Methods and Protocols*, ed. A. Navid (Totowa, NJ: Humana Press), 375–408. doi: 10.1007/978-1-61779-827- 6_13

Ploug, H. (2021). “SIMS and NanoSIMS techniques applied to studies of plankton productivity,” in *Research Methods of Environmental Physiology in Aquatic Sciences*, eds K. Gao, D. A. Hutchins, and J. Beardall (Singapore: Springer), 193–206. doi: 10.1007/978-981-15-5354-7_23

Polerecky, L., Masuda, T., Eichner, M., Rabouille, S., Vancová, M., Kienhuis, M. V. M., et al. (2021). Temporal patterns and intra- and intercellular variability in carbon and nitrogen assimilation by the unicellular cyanobacterium *Cyanothecae* sp. ATCC 51142. *Front. Microbiol.* 12:620915. doi: 10.3389/fmicb.2021.620915

Popa, R., Weber, P. K., Pett-Ridge, J., Finzi, J. A., Fallon, S. J., Hutcheon, I. D., et al. (2007). Carbon and nitrogen fixation and metabolite exchange in and between individual cells of *Anabaena oscillarioides*. *ISME J.* 1, 354–360. doi: 10.1038/ismej.2007.44

Schoffelen, N. J., Mohr, W., Ferdelman, T. G., Duerschlag, J., Littmann, S., Ploug, H., et al. (2019). Phosphate availability affects fixed nitrogen transfer from diazotrophs to their epibionts. *ISME J.* 13, 2701–2713. doi: 10.1038/s41396-019- 0453-5

Schoffelen, N. J., Mohr, W., Ferdelman, T. G., Littmann, S., Duerschlag, J., Zubkov, M. V., et al. (2018). Single-cell imaging of phosphorus uptake shows that key harmful algae rely on different phosphorus sources for growth. *Sci. Rep.* 8:17182. doi: 10.1038/s41598-018-35310-w
Schreiber, F., Littmann, S., Lavik, G., Escrig, S., Meibom, A., Kuypers, M. M. M., et al. (2016). Phenotypic heterogeneity driven by nutrient limitation promotes growth in fluctuating environments. *Nat. Microbiol.* 1:16055. doi: 10.1038/nnmicrobiol.2016.55

Stryhanyuk, H., Calabrese, F., Kümmel, S., Musat, F., Richnow, H. H., and Musat, N. (2018). Calculation of single cell assimilation rates from SIP-NanoSIMS-derived isotope ratios: a comprehensive approach. *Front. Microbiol.* 9:2342. doi: 10.3389/fmicb.2018.02342

Svedén, J. B., Adam, B., Walve, J., Nahar, N., Musat, N., Lavik, G., et al. (2015). High cell-specific rates of nitrogen and carbon fixation by the cyanobacterium *Aphanizomenon* sp. at low temperatures in the Baltic Sea. *FEMS Microbiol. Ecol.* 91:fiv131. doi: 10.1093/femsec/fiv131

Trembath-Reichert, E., Shah Walter, S. R., Ortiz, M. A. F., Carter, P. D., Girguis, P. R., and Huber, J. A. (2021). Multiple carbon incorporation strategies support microbial survival in cold subseafloor crustal fluids. *Sci. Adv.* 7:eabg0153. doi: 10.1126/sciadv.abg0153

Verity, P. G., Robertson, C. Y., Tronzo, C. R., Andrews, M. G., Nelson, J. R., and Sieracki, M. E. (1992). Relationships between cell volume and the carbon and nitrogen content of marine photosynthetic nanoplankton. *Limnol. Oceanogr.* 37, 1434–1446. doi: 10.4319/lo.1992.37.7.1434

Westfall, C. S., and Levin, P. A. (2017). Bacterial cell size: multifactorial and multifaceted. *Annu. Rev. Microbiol.* 71, 499–517.

Zimmermann, M., Escrig, S., Hübschmann, T., Kirf, M. K., Brand, A., Inglis, R. F., et al. (2015). Phenotypic heterogeneity in metabolic traits among single cells of a rare bacterial species in its natural environment quantified with a combination of flow cell sorting and NanoSIMS. *Front. Microbiol.* 6:243. doi: 10.3389/fmicb.2015.00243

**Conflict of Interest:** The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

**Publisher's Note:** All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Copyright © 2021 Polerecky, Eichner, Masuda, Zavřel, Rabouille, Campbell and Halsey. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.