To predict regimen efficacies for TB treatment and inform clinical trials, we consider the following complexities occurring during Mtb infection: spatial heterogeneities, Mtb phenotypic tolerance of drugs due to phenotype, drug-drug interactions, host-pathogen interactions, the host immune responses, and the effect of bacterial burden per granuloma (whereby density of bacteria affects PD, also known as inoculum effect in in vitro studies). In previous studies, we used our GranSim-PKPD model and tracked drug action using bactericidal assays performed in nutrient-rich conditions, macrophage assays or caseum mimics from various studies, together with INDIGO predictions for drug-drug interactions (Pienaar et al., 2015a; Pienaar et al., 2017; Budak et al., 2023; Pienaar et al., 2015b; Budak et al., 2024; Cicchese et al., 2021). Here, we improve our model of drug action by incorporating DiaMOND measurements performed under different in vitro growth conditions that model stressors that Mtb encounter in distinct granuloma microenvironments. These conditions have been shown to be predictive of in vivo and clinical outcomes, thereby improving our efficacy predictions (Larkins-Ford et al., 2021; Larkins-Ford et al., 2022). Further, we provide a methodology to calculate bactericidal activity, completely filling in the checkerboards, building on the “diagonal” of DiaMOND checkerboards, and single-drug responses. This methodology enables us to elucidate the picture of multi-scale drug-drug interactions and to assess optimal dosages to find combination regimens with the highest overall bactericidal activity, which, practically, can provide promising, mechanistically founded predictions to investigate in future studies. We name our new method GEODE: Granuloma Environment Optimization of Drug Efficacy (Figure 1).

Bacterial susceptibility to antibiotics can vary with bacterial burden, with higher bacterial burden potentially reducing drug efficacy (also called as “inoculum effect” in in vitro studies) (Brook, 1989), a phenomenon that has also been observed in Mycobacterium tuberculosis (Mtb) (Kirby and Dubos, 1947). To quantify and map the effect of Mtb burden, we conducted in vitro experiments to evaluate the bactericidal activity of antibiotics under a range of starting Mtb density. From these experiments, we generated pharmacodynamic (PD) data, modeled as Hill curves, which describe the drug-induced bacterial kill (measured as kill ratio in the GranSim framework; see Methods) as a function of antibiotic concentration for each bacterial load and drug tested (Figure 2).

We then translated this information into a pharmacokinetically relevant drug effect metric (E_PK), which was fitted to exponential decay curves, capturing how E_PK declines with increasing bacterial load (see Methods). As anticipated, incorporating an effect of varying levels of bacterial burden into the model resulted in a slower predicted rate of bacterial sterilization (Figure 3).

To validate our model capturing bacterial burden effect, we compare simulated sterilization curves to granuloma-scale data from marmosets (Budak et al., 2024) and macaques (Budak et al., 2023; Michael CT et al., 2025). We calculate sterilized granuloma ratios, the portion of a set of granulomas that sterilized, for both marmoset and macaque datasets, pooling granulomas from different animals and assuming a limit of detection as five colony forming units (CFU) per granuloma and, for some analyses, grouping granulomas by their pre-treatment CFU burden. We have used between 5 and 10 CFU as a threshold for our recent published works on granuloma modeling at all scales. This range is based on the calculation from NHP granulomas that are plated. This accounts for the amount of media that the granuloma is homogenized in. The model is flexible and the results are consistent within this range of CFU as guided by the NHP data (Lin et al., 2016). Additionally, consider that CFU is calculated as CFU/mL = (Number of Colonies Counted × Dilution Factor)/Volume Plated (in mL). So in this case, if we find 1 colony on the neat (most concentrated) plate: CFU/mL = (1 colony*1)/.2 mL = 5 CFU/mL => 5 CFU/granuloma in the ideal case. So both experimentally and theoretically five is a solid choice.

As marmosets tend to have progressive TB disease and also have many granulomas with uncontrolled CFU growth (Via et al., 2013) (see Figure 10B for GranSim vs marmoset CFU counts), we see that marmoset sterilized granuloma ratios best-match those from only high-CFU virtual granulomas (see Methods and Figure 10 for details of high- and low-CFU granulomas). We provide evidence that sterilized granuloma ratios from marmosets fit well to GranSim when we include representation of an effect of bacterial burden (Figures 3G–I). Further, sterilized granuloma ratios from macaques best match to low-CFU virtual granulomas, as most cynomolgus macaque granulomas tend to have controlled CFU growth (Maiello et al., 2018; Flynn et al., 2017). Here, we show that sterilized granuloma ratios from macaques agree with sterilization curves from low-CFU granulomas when we represent the effect of bacterial burden in our model (Figures 3D–F).

We developed a methodology to reconstruct the full checkerboard using the in vitro data measured via DiaMOND, enabling the estimation of bactericidal effects across all drug concentration combinations. This process involves three main steps: (1) converting growth-inhibition-based Hill parameters into bactericidal parameters; (2) deriving new Hill equation parameters for each drug combination using single-drug responses, equipotent combination data, and drug interaction scores (see Methods); and (3) calculating the killing rate for each combination based on these derived parameters. This approach provides a more comprehensive and scalable prediction of drug interactions than standard partial checkerboard methods (Figure 4).

As our goal in this work is to combine different drug screening approaches to improve regimen screening, in this work we combine in vitro combination response data for different bacterial growth conditions into GranSim for different drug regimens. In GranSim, we assume Mtb can assume three different phenotypes based on their distinct location: intracellular Mtb within macrophages, extracellular-replicating Mtb in granulomatous tissue and nonreplicating Mtb trapped within granuloma caseous-necrotic core. In GranSim, Mtb is modeled as existing in three distinct phenotypic states, reflecting well-characterized spatial niches within the granuloma: intracellular Mtb within macrophages, extracellular replicating Mtb in the surrounding tissue, and nonreplicating Mtb sequestered in the caseous-necrotic core. Using previously published in vitro data, we selected three in vitro models that together predicted treatment outcomes in the BALB/c relapsing mouse model (Larkins-Ford et al., 2021; Larkins-Ford et al., 2022): these included media with different sole carbon sources (butyrate or 0.2 mM cholesterol) and an acidic condition in standard growth medium. To capture drug pharmacodynamics (PD) relevant to specific Mtb subpopulations, we used the average killing rates from acidic and butyrate conditions to represent intracellular Mtb and the average rates from butyrate and high cholesterol to represent extracellular replicating Mtb. To model the PD of nonreplicating Mtb residing in caseous tissue, we developed a series of in vitro dormancy conditions varying in pH and oxygen levels: acidic normoxic (AN), neutral hypoxic (NH), and neutral normoxic (NN) (Methods and Table 1). Drug combinations were evaluated in these dormancy conditions using the DiaMOND methodology (Methods). We evaluated five clinically relevant combinations (HRZE, HRZM, RMZE, BPaL, and BPaMZ) under each of these dormancy conditions to evaluate their ability to predict clinical outcomes.

In vitro data measured in different growth stressors have been found to be a better predictor of in vivo outcomes than in vitro data from a single stressor alone (Larkins-Ford et al., 2021). Based on this prior work, we sought to generate composite stress condition profiles to parameterize GranSim. To do so, we first calculated the killing rate of each drug regimen under the individual stress conditions comprising the combination. We then aggregated these values using either the average, maximum or minimum killing rate, depending on the chosen strategy for modeling the combined effect, which results in the following combined conditions: AN_NH-min, AN_NH-max, AN_NH-mean, AN_NN-min, AN_NN-max, AN_NN-mean, NN_NH-min, NN_NH-max, NN_NH-mean, AN_NN-min, AN_NN-max, AN_NN_mean, AN_NN_NH-min, AN_NN_NH-max, AN_NN_NH-mean. For example, to aggregate all stress conditions in GranSim and to calculate the killing rate of a drug combination against Mtb considering the maximum of all conditions (AN_NN_NH-max), we determine the concentration of drugs an Mtb is exposed to on the simulation grid and calculate the killing rate of that drug combination considering each stress condition (AN, NN and NH) individually. We then take the maximum of killing rates in AN, NN and NH conditions to determine the killing rate of the aggregated condition AN_NN_NH-max (see Methods for details).

We first simulated clinically relevant regimens (HRZE, RMZE, HRZM, BPaL and BPaMZ) in GranSim considering all individual conditions and all possible aggregated conditions to determine the condition that is the best representative of clinical or in vivo outcomes (Figure 5). As the levels of nonreplicating bacteria trapped within caseum are likely the hardest population to clear with drugs due to drug tolerance and penetration issues, it is not surprising that PD of drugs acting on this population significantly affects sterilization dynamics. For example, BPaMZ has low efficacy assuming nonreplicating Mtb is exposed to certain growth conditions (AN_NH-min, NN_NH-min, AN_NN-min and AN_NN_NH-min, see Figures 5D,G,J,M), whereas BPaMZ sterilizes granulomas fastest when we assume nonreplicating Mtb is exposed to NN_NH-max, AN_NN-max and AN_NN_NH-max (Figures 5H,K,N). Moreover, the assumption that AN represents the microenvironment of nonreplicating Mtb results in a low sterilizing performance of HRZM (Figure 5A). By contrast, HRZM has one of the highest efficacies if one assumes that NN_NH-mean best describes that microenvironment (Figure 5I).

While sterilization dynamics of some regimens significantly depend on the chosen condition, some features of sterilization curves are common amongst most conditions. For example, RMZE performs the worst with almost all tested dormancy models except when BPaMZ has the lowest efficacy as in Figures 5D,G,J,M. In addition, BPaL is generally the best performer except when AN_NN-max is used for the PD of nonreplicating Mtb where BPaMZ sterilizes granulomas fastest (Figure 5K).

To determine growth conditions that reflect clinical and animal studies in a systematic way, we ranked the same five regimens for all simulated single- and composite dormancy conditions using Equation 24 (Figure 6). To choose which single- or composite dormancy condition we use to represent in vivo conditions, our primary criterion is that BDQ and PTM-containing regimens should outperform HRZE as the inclusion of BDQ and PTM has been previously shown to improve the efficacies of regimens compared to HRZE (Tasneen et al., 2016; Dawson et al., 2015; Cevik et al., 2024; Tasneen et al., 2011; Diacon et al., 2015). Only AN_NH-max, NN_NH-max, AN_NN-max, and AN_NN_NH-max fulfill these criteria, where BPaL and BPaMZ are the best-performing regimens (Figure 6). Next, we anticipate that HRZE should not outperform HRZM. The relative performances of HRZE and HRZM are debatable, as some studies claim HRZM has a higher efficacy than HRZE (Imperial et al., 2018; Nuermberger et al., 2004a; Burman et al., 2006; Conde et al., 2009), while clinical outcomes provide evidence that HRZM is not non-inferior to HRZE (Gillespie et al., 2014). Still, to our knowledge, there are no studies showing that HRZE could perform better than HRZM. These two criteria leave us with three conditions: AN_NH-max, NN_NH-max and AN_NN_NH-max. Since each of these conditions have similar sterilization dynamics, we continue our analyses by choosing one and mapping AN_NH-max to the PD of nonreplicating Mtb.

In the sections above, we isolate two new components critical to accurate antibiotic simulation - the effect of bacterial levels and a realistic representation of dormancy in GranSim. We now include a third feature: in vitro measurements of drug-drug interactions (see Methods) to complete our expanded drug simulation platform: GEODE (Granuloma Environment Optimization of Drug Efficacy). GEODE is our novel, drug regimen efficacy prediction model that combines drug-drug interaction PD with our GranSim computational model of PK/PD dynamics within a granuloma. To validate GEODE, we compare our treatment simulations to both NHP and clinical studies. We use datasets from two different types of NHPs: macaques and marmosets. Cynomolgus macaque are an accurate model of human Mtb infection, whereas common marmosets almost exclusively develop progressive disease (Via et al., 2013; Flynn et al., 2017). We analyzed three features of simulations: (a) temporal CFU trajectories within all granulomas in response to treatment, (b) rankings of regimens based on granulomas sterilization times, and (c) change in overall granuloma metabolic activity after treatment, which can be simulated and has been experimentally quantified via FDG avidity from PET/CT scans (Budak et al., 2023; Lin et al., 2013; White et al., 2017). For (a–c), data from the NHPs are only available for three of our five clinically relevant treatment regimens, namely, HRZE, RMZE and BPaL. Therefore, our comparisons are only for that.

GranSim is a computational, hybrid, multi-scale, agent-based model that we have developed and continuously curated with NHP data for 2 decades; GranSim simulates host-pathogen interactions upon Mtb infection in an area of lung tissue (6 mm × 6 mm) (Segovia-Juarez et al., 2004; Ray et al., 2009; Fallahi-Sichani et al., 2011; Cilfone et al., 2013). GranSim contains cellular agents, such as immune cells (e.g., macrophages, T cells), Mtb and signaling molecules (chemokines and cytokines). We define the rules that describe the cellular and molecular processes and the agent interactions within GranSim based on experimental observations from NHP studies (for a complete description of our rules, see http://malthus.micro.med.umich.edu/GranSim/). We represent these rules using mathematical equations containing parameters and we calibrated all parameters to NHP datasets (Segovia-Juarez et al., 2004; Ray et al., 2009; Cilfone et al., 2013). Each simulation begins with a single infected macrophage in the center of the lung tissue grid that triggers recruitment of immune cells to the infection site. In simulations, structured granulomas are a consistent emergent result of the Mtb and immune dynamics, and that structure shares similar spatial features as NHP and human granulomas—e.g., a necrotic core, a cellular rim composed of macrophages surrounding the necrotic core, and a T-cell cuff as the outermost ring.

To capture the heterogeneity of Mtb metabolic states, we included three types of Mtb based on their location within a granuloma: intracellular Mtb reside within macrophages, nonreplicating Mtb are those trapped in the necrotic core, and replicating extracellular Mtb are those within granulomatous tissue (i.e., the non-necrotic region). These Mtb phenotypes have distinct features in terms of replication rates, interactions with immune cells and susceptibility to anti-TB drugs (Lenaerts et al., 2015).

To analyze regimen efficacies at the granuloma scale, we use GranSim to generate an in silico granuloma library that contains 200 granulomas. We generate this granuloma library by sampling 250 parameter sets within biologically relevant ranges using Latin Hypercube Sampling and simulate granulomas parameterized by each parameter sample for 300 days post-infection including three replications of each (to address aleatory (stochastic) uncertainty (Marino et al., 2008), 750 simulations in total. The parameter ranges are calibrated to CFU counts from macaque granulomas (Figure 10A). We then categorize granulomas into low- and high-CFU granulomas (as described in the next paragraph) and randomly pick 100 low- and 100 high-CFU granulomas from the set of 750 granulomas, having 200 granulomas in total. Due to a random selection process, it is possible to have granulomas with the same parameter set but different stochastic processes (i.e., different seeds that determine stochastic events) within the in silico granuloma library.

To assess how effective regimens are in treating different types of granulomas, we grouped them by their capacity to control Mtb growth. To do this, we grouped virtual granulomas into two groups based on their bacterial growth trends after day 300 post-infection: low-CFU and high-CFU granulomas. LTBI patients mainly have sterile or low-CFU granulomas with stable CFU levels, whereas active TB patients likely have at least one high-CFU granuloma with uncontrolled CFU growth (Cicchese et al., 2020a; Lin and Flynn, 2018; Gideon et al., 2015). As in our previous work (Budak et al., 2023; Budak et al., 2024), we define granulomas as low-CFU if they have nonzero and less than 10⁴ CFU at the end of 300 days post-infection, and their CFU levels do not increase more than 50 CFU in the last 20 days of simulation (blue curves in Figure 10). We define granulomas as if they have CFUs between 10⁴ and 10⁷ at the end of the 300-day simulation or their CFUs increase more than 50 CFU in the last 20 days of simulation (red curves in Figure 10) (Our model is flexible and can easily update these ranges). Out of 750 simulations, we randomly selected 100 low-CFU and 100 high-CFU granulomas, generating a library of 200 granulomas. Note that marmosets are very small animals and GranSim timing was based mostly on macaques and human datasets; thus, infection-dynamic timing in marmosets is different than simulations. In fact, CFU counts from marmoset granulomas correlate better with later time points in infection dynamics as simulated by GranSim (Figure 10B).

We previously included a PK model within GranSim that simulates spatial distributions of antibiotics within simulated granulomas following administration of antibiotics (Pienaar et al., 2015a; Pienaar et al., 2017; Pienaar et al., 2015b), which we continue to use in the present work. Briefly, we first simulate plasma PK using a compartmental model with a system of ordinary differential equations. Drugs concentrations in plasma permeate into lungs, which we model via flux terms of drugs from plasma onto the simulation lung grid through vascular source agents. Once a drug permeates into simulated lung tissue, it diffuses through cellular tissue, binds to macromolecules (e.g., caseum, epithelial tissue, etc.), and partitions into macrophages (see http://malthus.micro.med.umich.edu/GranSim/ for the full list of equations and code files). We calibrate the parameters that define PK processes using in vivo studies from humans or animal models (Cicchese et al., 2020a; Budak et al., 2023; Budak et al., 2024). GranSim currently represents eight antibiotic drugs used for treatment of TB: isoniazid (H or INH), rifampicin (R or RIF), pyrazinamide (Z or PZA), ethambutol (E or EMB), moxifloxacin (M or MXF), Bedaquiline (B or BDQ), Pretomanid (Pa or PTM), Linezolid (L or LZD).

In vitro assays were performed using DiaMOND experimental design in Mycobacterium tuberculosis (Erdman strain). In vitro data for the cholesterol, butyrate, and acidic conditions were obtained from Larkins-Ford et al. (2021). For the dormancy (non-replicating) conditions in this study (Table 1), M. tuberculosis was initially cultured in Middlebrook 7H9 broth supplemented with OADC and Tween-80, then acclimated to a lipid-rich medium containing oleic acid, palmitic acid, arachidonic acid, linoleic acid, stearic acid, and cholesterol.

The dormancy media were buffered to a pH of 5.5 for acidic conditions and pH 7.0 for neutral conditions using MES (2-(N-morpholino)ethanesulfonic acid) and MOPS (3-(N-morpholino)propanesulfonic acid), respectively. Normoxic cultures were maintained under ambient oxygen (vented flasks), while hypoxic cultures were incubated in low-oxygen conditions (non-vented flasks). The dormancy cultures were incubated at 37 °C for 4 weeks prior to drug treatment, followed by 1 week of drug exposure. Given the low metabolic state of non-replicating cells, cultures were plated on charcoal agar (Gold et al., 2016) to inactivate residual antibiotics and promote outgrowth. Luminescence was measured to assess bacterial survival, as the Mtb strain used in this study carries the autoluminescent reporter plasmid pMV306hsp+LuxG13 (Larkins-Ford et al., 2021). All assays were performed in biological triplicate, and the median value across replicates was used for downstream analyses (see Supplementary File 2 for the raw data).

We create a methodology for predicting bactericidal activities of various drug dose combinations. Note that we distinguish between two major measurements of the potency of a drug: bactericidal activity and bacterial inhibition. To infer the bactericidal activity of multi-drug combinations in GEODE, we use data provided from in vitro assays (Cokol et al., 2017; Larkins-Ford et al., 2021) in the form of fractional inhibitory concentration values (FIC, see Glossary). Concretely, inhibitory results are represented by a Hill function that defines the level of bacterial growth inhibition as a function of drug concentration. The general form of a Hill curve can be represented as E C = Emax C h C h + C 50 h (1) where E is the effect of a drug that depends on its concentration C , Emax is the maximum effect a drug can reach, C 50 is the concentration needed to reach the half-maximum effect ( Emax /2) and h is the Hill constant, which characterizes the steepness at the inflection point. Here, we present two sets of Hill curve parameters: growth inhibition-based parameters ( Emax i n h ; C 50 inh , h inh ) obtained from in vitro assays and, inferred from that in the next subsection, bactericidal activity-based parameters ( Emax frac _ killed , i , Emax i , C 50 N , i , h i for drug i ; Emax c o m b and Emax frac _ killed , comb for combinations).

Datasets from these assays represent growth inhibition-based PD for single-drug dynamics (x- and y-axis of a standard checkerboard assays) and equipotent drug combinations (diagonals of the same checkerboard) (Larkins-Ford et al., 2021). Herein we derive a method to calculate the bactericidal activities of all possible dose combinations (i.e., the whole checkerboard) by using growth-inhibition-based dose response curves for single drugs and the equipotent combinations, and colony counts for different growth inhibition levels.

Using GranSim, we simulated each regimen using our in silico granuloma library containing 100 high-CFU and 100 low-CFU granulomas, totaling 200 granulomas, for a maximum of 180 days (shorter simulations to be consistent when comparing to experimental/clinical studies are clarified below). Unless otherwise stated, we administer human or human-equivalent doses of each drug within a regimen daily, based on the model we used for PK calibration (Budak et al., 2023; Budak et al., 2024): 5 mg/kg of isoniazid (World Health Organization, 2014), 10 mg/kg of rifampicin (World Health Organization, 2014), 7 mg/kg of moxifloxacin (World Health Organization, 2014), 20 mg/kg of ethambutol (World Health Organization, 2014), 25 mg/kg of pyrazinamide (World Health Organization, 2014), 20 mg/kg of Bedaquiline (human-equivalent dose for rabbits), 20 mg/kg of pretomanid (human-equivalent dose for rabbits), and 90 mg/kg of linezolid (human-equivalent dose for rabbits) (Budak et al., 2024).

We performed three complete sets of antibiotic treatment simulations:

a. Treatment simulations to replicate clinically relevant regimens and to choose the right stress condition that models the PD of nonreplicating Mtb:We simulated regimens that have been widely studied in clinical and preclinical trials, namely, BPaMZ (Cevik et al., 2024; Tweed et al., 2019; Li et al., 2017; Xu et al., 2019), BPaL (Tasneen et al., 2016; Xu et al., 2019; Bigelow et al., 2020; Berry et al., 2022; Conradie et al., 2020; Conradie et al., 2022), HRZE (Burman et al., 2006; Gillespie et al., 2014; Dorman et al., 2009; Johnson et al., 2009), RMZE (Nuermberger et al., 2004a; Gillespie et al., 2014; Nuermberger et al., 2004b; Li et al., 2015; Dorman et al., 2009) and HRZM (Nuermberger et al., 2004a; Burman et al., 2006; Conde et al., 2009; Gillespie et al., 2014; Li et al., 2015), for 180 days using 15 different sets of in vitro measurements to model the PD and drug interaction against nonreplicating Mtb (Table 1).

A regimen set is a group of regimens to be ranked based on their simulated efficacies using GranSim. We assess regimen efficacies based on sterilization times of each granuloma and rank them accordingly. To rank the regimens in a regimen set, we calculate a ranking score for each regimen as described in Budak et al. (2024). We refer to the reference regimen as the regimen for which the ranking score is calculated. Briefly, we compare the sterilization times of the reference regimen pairwise to that of the rest of the regimens in the regimen set. We do this using a one-tailed t-test, assuming the measures are normally distributed. A regimen is considered more potent than a reference regimen if its sterilization times are significantly lower than that of the reference regimen, and vice versa (p < 0.05, one-tailed t-test). We then count the number of significantly more (and less) potent regimens than the reference regimen to calculate ranking score as: Ranking   score reference   regimen = no .  of more potent  regimens than the  reference regimen − no .  of less  potent regimens than  the reference regimen (24)

We calculate ranking scores of each regimen in a regimen set and used these scores to rank regimens, in that higher ranking scores and lower ranks mean the regimen performs better overall within that specific regimen set.

To assess the impact of bacterial burden on drug killing rates (i.e., inoculum effect in in vitro studies), M. tuberculosis was cultured under a single-stressor assay, the butyrate condition. Cultures were grown to mid-log phase, then diluted to three starting bacterial load levels: OD₆₀₀ of 0.05, 0.1, and 0.15 for the drug susceptibility assay (see Supplementary File 3 for the raw data). From these experiments, we derive Hill curve parameters ( E max , C 50 and h ) that yield drug kill rates dependent on drug concentrations for varying drugs (d) and bacterial load (OD600) (Figure 2). We define the drug effect as kill ratio, which indicates the ratio of the killed bacteria to total bacteria (i.e., the portion of the killed bacteria) depending on the drug concentration: Kill ratio C = E max C h C h + C 50 h (25)

Since the PK of each drug varies and the drug effect depends on drug concentration, we define PK-relevant drug effect E PK for each bacterial load and drug as: E PK OD 600 , d = Emax OD 600 , d ∫ 0 Cmax d C h OD 600 , d C h OD 600 , d + C 50 OD 600 , d h OD 600 , d dC Cmax d (26) where E ma x OD 600 , d , C 50 OD 600 , d and h OD 600 , d are Hill curve parameters for drug d and bacterial load value of OD600, C ma x d is the maximum concentration of drug d based on GranSim PK simulations (Supplementary Table S11; Supplementary Figure S4, see Figure 2 for the Hill curves of all drugs and bacterial loads).

We then convert OD600 values from bacterial culture to CFU concentrations (to make it consistent with GranSim which is based on CFU): CFU / ml OD 600 = 1.198 · 10 8 / ml · OD 600 + 1.297 · 10 6 / ml (27)

OD-based CFU quantification (Equation 27) was based on standard microbiology technique where replicate Erdman-strain Mtb were grown at 37 °C/5% CO₂ in liquid culture without agitation. The OD600 was measured with a spectrophotometer at different timepoints, and then dilutions were plated on 7H11 plates so that the CFU/mL could be calculated and compared against the OD from that culture.

Once we calculate CFU/mL values and their corresponding E PK values for each drug using Equation 26, we then fit these datapoints using an exponential function in the form of A e − B x such that: E PK CFU ml , d = A d e − B d · CFU ml (28) where A d and B d are drug-dependent parameters (see Table 2; Figure 12).

In GranSim, we calculate a modified killing rate k bacterial   l o a d by scaling the killing rate k (calculated in 5) using E P K values depending on CFU concentration: k bacterial   load = k · E PK CFU ml G r a n S i m E PK CFU ml OD 600 = 0.05 (29) where CFU ml G r a n S i m is the CFU concentration of simulated granulomas. In Equation 29, we normalize it to E PK CFU ml OD 600 = 0.05 , as in vitro experiments (that were used to calculate k ) were performed using Mtb with OD 600 = 0.05 .

When multiple drugs are available in a microgrid, g , in GranSim, then we calculate the bacterial load-dependent killing rate k at grid point g by taking the mean of k bacterial   l o a d values of all available drugs: k bacterial   load g = k · ∑ d = 1 n A d e − B d · CFU ml GranSim n · E PK CFU ml OD 600 = 0.05 (30) where n is the number of drugs available in the microgrid g .

Adult male Mauritian cynomolgus macaques were used in this study and obtained from Bioculture US. All experimental procedures involving care of animals complied with ethical regulations of the University of Pittsburgh School of Medicine Institutional Animal Care and Use Committee. Macaques were housed and cared for in accordance with local, state, federal, and institute policies in facilities accredited by the American Association for Accreditation of Laboratory Animal Care, under standards established in the Animal Welfare Act and the Guide for the Care and Use of Laboratory Animals as mandated by the U.S. Public Health Service Policy. Macaques were monitored for physical health, food consumption, body weight, temperature, complete blood counts, and serum chemistries. Examination of animals was performed in quarantine upon arrival at the University of Pittsburgh to assess physical health and confirm no previous Mtb infection was detected via ELISpot assays. All Mtb infections were performed a Biosafety Level 3 (BSL3) facility. Bronchoscopic instillation with Mtb Erdman was performed as previously described with an average dose of ∼19 CFU (Capuano et al., 2003; Lin et al., 2009). Veterinary staff regularly monitored clinical signs following challenge, including appetite, behavior and activity, weight, erythrocyte sedimentation rate, Mtb growth from gastric aspirate, and coughing.

Five adult male Mauritian cynomolgus macaques were obtained from Bioculture, Inc (Mauritius), with an age range of 5–9 years (Supplementary Material 1 provides all macaque data from the BPaL study). Macaques were co-housed in a Biosafety Level 3 facility and infected with 8–22 CFU M. tuberculosis strain Erdman via bronchoscope as described previously (Winchell et al., 2023). Approximately 8 weeks after M. tuberculosis infection, the animals were treated with the BPaL regimen daily via oral administration in food treats for 4 weeks. Drug doses: Bedaquiline (B) 40 mg/kg, Pretomanid (P) 62 mg/kg, Linezolid (L) 20 mg/kg. RMZE and HRZE treated animals were obtained from Valley Biosystems (5 adult males with age range of 4–9 years) and were infected with 8–21 CFU (Budak et al., 2023). Approximately 12 weeks after Mtb infection, the animals were treated with HRZE (n = 3) or RMZE (n = 2) for 8 weeks. Drug doses: Isoniazid (H) 15 mg/kg; Rifampicin (R) 20 mg/kg; Pyrazinamide (Z) 150 mg/kg; Ethambutol (E) 55 mg/kg; Moxifloxacin (M) 35 mg/kg. For all drug regimens, PET CT scans were performed immediately prior to initiation of drug treatment, and 4 and 8 weeks after initiation of drug treatment using 18-F fluorodeoxyglucose as a PET probe as previously described (White et al., 2017; Winchell et al., 2020). The BPaL-treated animals had only the 4 weeks post-drug treatment scan directly prior to necropsy with the shorter treatment time. The HRZE- and RMZE-treated animals had both a 4 weeks and 8 weeks scans. Granuloma size (measured by maximum diameter in mm) and avidity (measured by maximum SUV and corrected for granuloma size) were evaluated using OsiriX DICOM imaging software in the axial view. Avidity was measured for granulomas resulting in SUVR values pre- and post-drug treatment (White et al., 2017). A fold change was calculated per granuloma to measure drug effect on an individual granuloma basis (see Equation 31 below). At necropsy the final PET CT scan was used as a map to identify and individually isolate all granulomas or other pathologies, uninvolved lung lobes, thoracic lymph nodes and any signs of extrapulmonary disease. Each sample was plated for bacterial burden individually on 7H10 plates; the plates were incubated at 37° with 5% CO2 and counted for colony forming units (CFU) after 3 weeks. The fraction (%) of sterile (CFU negative) granulomas was calculated for each animal. Fold change in SUVR = log 2 prenecropsy   SUVR + 0.01   predrug   SUVR + 0.01 (31)