Mechanistic Modeling of cfDNA Fragmentome
Dynamics Predicts Progression to Immunotherapy

Linh Nguyen Phuong, Frédéric Fina, Laurent Greillier, Pascale Tomasini,
Jean-Laurent Deville, Audrey Boutonnet, Frédéric Ginot, Jean-Charles
Garcia, Sebastien Salas, Sébastien Benzekry

► To cite this version:

Linh Nguyen Phuong, Frédéric Fina, Laurent Greillier, Pascale Tomasini, Jean-Laurent Deville, et al.  
Mechanistic Modeling of cfDNA Fragmentome Dynamics Predicts Progression to Immunotherapy.  
2025. hal-05241421

HAL Id: hal-05241421

https://inria.hal.science/hal-05241421v1

Preprint submitted on 20 Jan 2026

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished 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_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
1 Mechanistic Modeling of cfDNA Fragmentome Dynamics Predicts

2 Progression to Immunotherapy

3 Linh Nguyen Phuong¹,², Frédéric Fina³, Laurent Greillier¹,²,⁴,⁵, Pascale Tomasini⁴,⁵,

4 Jean-Laurent Deville⁵, Audrey Boutonnet⁶, Frédéric Ginot⁶, Jean-Charles Garcia⁶,

5 Sebastien Salas¹,²,⁴,* , Sebastien Benzekry¹,²,*

6

7 ¹COMPutational pharmacology and clinical Oncology, Centre Inria d'Université Côte

8 d'Azur, Sophia-Antipolis, France

9 ² Cancer Research Center of Marseille, Institut Paoli-Calmettes, Inserm UMR1068,

10 CNRS UMR7258, Aix Marseille University UM105, Marseille, France

11 ³ID-Solutions oncology, Marseille, France

12 ⁴Multidisciplinary Oncology & Therapeutic Innovations Department, Assistance

13 Publique—Hôpitaux de Marseille,

14 ⁵Assistance Publique—Hôpitaux de Marseille, Timone Hospital, Aix Marseille

15 University, Marseille, France

16 ⁶Adelis Technologies, Labège, France

17

18 * = joint senior authors
19 ABSTRACT

20 Background

21 Plasma cell-free DNA (cfDNA) shows promise as predictive cancer biomarker, but
22 mechanisms governing cfDNA production/fragmentation/elimination dynamics, and
23 their relationships with tumor burden and disease progression, remain poorly
24 understood.

25 Methods

26 We developed a mechanistic model jointly describing short (75-<580 bp) and long
27 (≥580-1650 bp) cfDNA dynamics alongside tumor kinetics in 112 advanced cancer
28 patients receiving immune checkpoint inhibition, using a population approach.

29 Results

30 The model successfully described complex cfDNA patterns, including treatment-
31 initiation spikes. It revealed large inter-patient variability in kinetic cfDNA parameters,
32 and a 7.4-fold higher shedding rate for short versus long fragments. A model-derived
33 parameter from 6-weeks data — reflecting enhanced release or reduced elimination of
34 short fragments — was significantly associated with progression-free survival (PFS)
35 (HR=1.6 [95% confidence interval (CI): 1.2-2.2], p = 0.001). Notably, adding this
36 parameter to baseline clinical prognostic variables improved PFS prediction (C-index
37 0.78 [95% CI: 0.73-0.89] vs 0.80 [95% CI: 0.74-0.90], p < 0.0001).

38 Conclusion

39 Our model provides quantitative insights into cfDNA biology and offers a non-invasive
40 way to monitor and predict resistance before imaging.

41

42 **Keywords:** cell-free DNA; fragmentome; size profile; biomarker; mechanistic
43 modeling; immune-checkpoint inhibition.

44
45 Author Summary

46 Using an ordinary differential equation system, we modeled the kinetics of both short
47 (75–580 bp) and long (580–1650 bp) cfDNA fragments during immunotherapy,
48 assuming that variations in cfDNA largely arise from tumor-derived release through
49 chaotic apoptosis, necrosis, and active secretion. Using a nonlinear mixed-effect
50 statistical framework on clinical data, this first mechanistic model of the joint size-
51 dependent cfDNA and tumor dynamics was identifiable, provided quantitative insights
52 of both short and long fragment biology, and captured complex cfDNA patterns,
53 including the initial spike often observed at treatment onset, which may reflect rapid
54 immunotherapy-induced tumor shrinkage.

55 We then applied a landmark analysis strategy using data available at six weeks to
56 assess associations with progression-free survival (PFS). Notably, we identified a
57 significant association between the short-fragment release rate and shorter PFS.
58 Overall, this mechanistic modeling approach shows promise for leveraging longitudinal
59 cfDNA fragment sizes to monitor immunotherapy response, particularly to guide early
60 treatment-discontinuation decisions.
61 INTRODUCTION

62 Through the last decades, circulating cell-free DNA (cfDNA) has demonstrated high
63 potential as a promising biomarker for treatment monitoring through liquid biopsy
64 approaches. This methodology offers several advantages over traditional tissue-based
65 assessments: it provides access to systemic biological information, catching both
66 tumor and wild-type data, reflects intra- and inter-tumor genetic heterogeneity, and
67 enables frequent, non-invasiveness, and cost-effective sampling throughout treatment.
68 Study of cfDNA fragmentation patterns, known as fragmentomics¹, has emerged as an
69 important area of research for understanding tumor biology and immune system
70 function, while offering clinical applications for disease detection and therapy
71 monitoring². In particular, the size distribution profile of cfDNA fragments presents
72 interesting characteristics, with sharp peaks at the base pair positions corresponding
73 to k-nucleosomes (e.g., released fragments from mono-nucleosomes have a typical
74 size of 166 base pairs)³.

75 In oncology, immune checkpoint inhibitors (ICI) have revolutionized cancer treatment,
76 significantly improving survival outcomes for patients with advanced malignancies over
77 the last two decades. However, 20 to 40% of patients fail to respond to
78 immunotherapy⁴, experiencing disease progress, or developing immune-related
79 adverse events. This therapeutic variability has made the understanding of biological
80 dynamics of ICI response and its prediction early on-treatment the most pressing
81 challenges in the last years.

82 Conventional approaches to monitor ICI response face several critical limitations.
83 Current imaging-based assessments through iRECIST criteria (Response Evaluation
84 Criteria in Solid Tumors)⁵ provide static evaluations, preventing quick individualized
85 therapeutic adaptation. Established predictive biomarkers, including Programmed
86 Death-Ligand 1 (PD-L1) protein expression and Tumor Mutation Burden (TMB),
87 require invasive tissue biopsies, that cannot be performed frequently. Furthermore,
88 these methods demonstrate limited standardization and poor predictive accuracy of
89 resistance, with positive predictive values (PPV) of 34 and 42%, respectively, in multi-
90 cancer contexts⁶. In contrast, the cfDNA fragmentome offers a non-invasive, cheap
91 and promising source of biomarkers to predict and monitor response / resistance.
Most investigations have concentrated on pre-treatment cfDNA analysis to enable early treatment adjustments⁷. While many studies have focused on circulating tumor DNA (ctDNA) characteristics, such as blood TMB, copy number instability, mutation burden⁸⁻¹⁰, several notable studies have accounted for global cfDNA parameters, such as global cfDNA concentration¹¹, or fragmentome data¹²,¹³. Recently, Janke et al. found association between baseline fragment length alterations and PFS, outperforming and complementing PD-L1 protein expression in advanced non-small cell lung cancer (NSCLC) patients under immunotherapy¹⁴. Several studies have also assessed overall survival (OS) outcomes based on cfDNA features collected during treatment, including both absolute and relative to baseline measurements¹⁵⁻¹⁸. Notably, Assaf et al. were among the first to demonstrate a strong association between early changes of ctDNA as well as total cfDNA concentration and OS (7.53 vs. 25.2 months in the progressors vs. stable/responder groups, C-index = 0.66) in NSCLC patients, using a machine learning framework¹⁹.

While previous studies have highlighted pre- and on-treatment cfDNA markers, they treated these as static measurements rather than modelling the continuous temporal dynamics throughout treatment. Mechanistic modeling offers a powerful approach to represent underlying biological processes over time, enabling not only physiological interpretation and assessment of biological assumptions, but also providing computational markers potentially predictive of outcome.

The first mechanistic models describing systemic quantities beyond pharmacokinetics in oncology were tumor growth inhibition models²⁰. These relied on the successful population approach widely employed in pharmacometrics²¹. Modeling of other circulating biomarkers has also been proposed, including for treatment prediction ²²,²³. For instance, Netterberg et al. linked mean circulating tumor cell (CTC) counts to tumor growth in patients with metastatic colorectal cancer²⁴, showing that CTC predictions from tumor size measurements via Poisson model outperformed tumor kinetics (TK) as OS predictors. In a recent study, we demonstrated that modeling kinetics of simple laboratory blood markers could enhance prediction of OS compared to TK alone²⁵. However, few empirically-validated models have yet been proposed for cfDNA data. Avanzini et al. first employed mechanistic modeling to determine an optimal blood-based screening frequency for cancer diagnosis based on cfDNA data ²⁶. Janssen et al. developed a nonlinear mixed-effects model to describe EGFR mutation
concentrations in NSCLC patients receiving chemotherapy. The model's predicted concentrations were found to be predictive of disease progression during treatment²⁷.

Li et al. recently proposed a stochastic approach modeling ctDNA release by tumor cells and elimination under targeted therapy or radiotherapy in NSCLC patients²⁸.

Ribba et al. developed the first joint ctDNA—TK model to assess atezolizumab response in NSCLC patients, incorporating a coefficient to link tumor size decay rate to the ctDNA growth rate²⁹. Beyond the mere ctDNA, limited research has focused on mechanistic modeling of the total cfDNA, and even less on fragment size-structured dynamics. Using Markov chain modeling to describe cfDNA fragmentation and simulate size distributions, Tsui et al. demonstrated that the elimination of cfDNA fragments depends on their size, with smaller fragments being cleared more easily from the circulation than longer ones³⁰. However, no existing model has integrated size-dependent cfDNA concentrations with TK.

This study presents the first mechanistic model of joint TK – size-structured cfDNA dynamics. We relied on longitudinal cfDNA size distribution data by grouping the fragments into two groups: short fragments – corresponding to mono-, di-, and tri-nucleosomal fragments and mainly resulting from apoptosis³ – and long fragments mainly resulting from necrosis³¹. We had two primary objectives: 1) to quantitatively elucidate the biological mechanisms underlying cfDNA release and elimination and 2) to assess the predictive ability for progression of early, on-treatment data of the dynamics of the cfDNA fragmentome.
146 **METHODS**

147 **Study design and patient population**

148 The Size CfDNA Immunotherapy Signature Monitoring (SChISM) clinical study has
149 been introduced and comprehensively described previously³². Briefly, SChISM is a
150 prospective, non-interventional clinical study that aims to investigate cfDNA as a
151 predictive tool in 337 patients with metastatic and/or recurrent NSCLC, head and neck
152 squamous cell carcinoma (HNSCC), clear-cell renal cell carcinoma (ccRCC), urothelial
153 carcinoma (UC), and melanoma initiating treatment with ICIs as standard care. cfDNA
154 samples were collected before each cure of ICI until the first imaging evaluation, then
155 every two cycles for one year, and at recurrence. The study cohort was limited to
156 patients with complete data and adequate follow-up for longitudinal modeling.
157 Exclusions included melanoma diagnosis (i.e., non-carcinoma, n = 110), non-
158 evaluable disease status (n = 43), insufficient follow-up data (n = 58), missing tumor
159 size data (n = 10), and elevated cfDNA concentration (> 1000 pg/µL at one timepoint)
160 (n = 4).

161 **Blood collection and cfDNA quantification**

162 Blood samples (Roche cfDNA Collection Tube) were processed through two
163 centrifugation steps to isolate plasma and remove cellular debris. The plasma was then
164 aliquoted and frozen at -80°C for storage. Plasma was analyzed by Adelis
165 Technologies (Labège, France), using the BIABooster device for size profile
166 characterization of cfDNA. cfDNA was analyzed using capillary electrophoresis with
167 fluorescence detection. The system concentrated DNA fragments while removing salts
168 and proteins, then separated fragments by size (75-1650 base pairs (bp) range) with
169 10 fg/µL sensitivity³³.

170 **Variables**

171 Clinical and biological data were prospectively collected for each patient included in
172 the study. Collected variables included age, sex, Eastern Cooperative Oncology Group
173 (ECOG) performance status, tumor type (NSCLC, HNSCC, ccRCC, or UC), and
174 disease stage at the time of immunotherapy initiation (metastatic and/or recurrent).
175 Treatment response was assessed every three months using contrast-enhanced
176 computed tomography (CT) scans, interpreted according to iRECIST criteria⁵. The sum
177 of the largest diameters (SLD) of target lesions was computed from these CT scans.
178 The date of radiological progression, when applicable, was recorded, as well as the
179 date of last follow-up or death for survival analyses. In addition, routine laboratory
180 parameters were collected at baseline and during follow-up. The neutrophil-to-
181 lymphocyte ratio (NLR) was calculated for each patient and selected as a nonspecific
182 serum biomarker that has been extensively studied as both a prognostic and predictive
183 factor in the context of immunotherapy³⁴,³⁵.

184 One cfDNA curve is defined as the concentration of fragments in pg/µL according to
185 the fragment size in bp (Figure 1A). Three quantitative variables were derived from
186 each curve for longitudinal modeling, based on different ranges of size:

187 - The total concentration *D*, defined as the area under the cfDNA curve between 75
188 and 1650 bp and computed by the trapezoid method.

189 - The concentration in pg/µL of short cfDNA fragments measuring between 75 and
190 580 bp, denoted as *D*<sub>S</sub>.

191 - The concentration in pg/µL of short cfDNA fragments measuring between 580 and
192 1650 bp, denoted as *D*<sub>L</sub>.

193 The threshold between short and long fragments was arbitrarily set at 580 bp. It
194 corresponds to the boundary after the first three peaks in cfDNA size profiles, which
195 represent DNA fragments typically released through apoptosis³.
196

Moreover, this threshold allowed to distinguish between two quantities that vary
197 differently over time, enabling meaningful analysis of their distinct patterns.

198 Longitudinal modeling

199 A mechanistic model was developed to jointly describe the dynamics of global cfDNA
200 concentration short fragments Ds and long fragments Dl of cfDNA with tumor size T under
201 immunotherapy (Figure 1B). Short fragments arise from apoptosis³, while longer
202 fragments (genome-wide) originate from necrosis³¹, though these are less abundant in
203 plasma³. CfDNA clearance occurs primarily through hepatic, splenic, and renal
204 mechanisms³⁶,³⁷. However, the precise elimination kinetics remain poorly understood,
205 with conflict evidence regarding cfDNA clearance saturation³⁶ and its dependence with
206 fragment sizes³⁰. Active secretion processes can release fragments of variable sizes
207 (~100 bp to >10 Kbp)³⁸ via extracellular vesicle formation for cellular homeostasis
208 maintenance, or through micronuclei release during replication error correction,
209 particularly in patients with chromosomal instability³⁹.
210 Thus, the model relied on the following biological assumptions:

211 - tumor cells (assimilated to the sum $T(t)$ of largest diameters of all lesions)
212 comprise two sub-populations: cells sensitive to treatment ($T_S$) and resistant ones
213 ($T_R$), governed by first ordered kinetics driven by a shrinkage parameter $\beta$ and a
214 regrowth parameter $\alpha$.

215 - Short fragments ($D_S$) are proportionally released, via a rate $\lambda_S$ (expressed in pg.µL⁻¹.mm⁻¹), through:

216 o Cancer cells ($T$) proliferation, through active secretion.

218 o Treatment-sensitive cancer cells ($T_S$) death, through apoptosis³.

219 - Long fragments $D_L$ are proportionally released, via a rate $\lambda_L$ (expressed in pg.µL⁻¹.mm⁻¹), through:

221 o Cancer cells ($T$) proliferation, due to active secretion and/or necrosis of the tumor microenvironment.

223 o Treatment-sensitive cancer cells ($T_S$) death, through necrosis.

224 - Short fragments are linearly cleared from the circulation at a clearance rate $k_{D_S}$ (expressed in month⁻¹).

226 - Long fragments are linearly cleared from the circulation at a clearance rate $k_{D_L}$ (expressed in month⁻¹).

$$ \frac{dT_R}{dt} = \alpha \cdot T_R $$ 
$$ \frac{dT_S}{dt} = \begin{cases} \alpha \cdot T_S, & \text{if } t < 0 \\ (\alpha - \beta) \cdot T_S, & \text{if } t \ge 0 \end{cases} $$ 
$$ T_{R0} = T_R(t_0) $$ 
$$ T_{S0} = T_S(t_0) $$ 
$$ T = T_S + T_R $$

$$ \frac{dD_S}{dt} = \lambda_S \cdot (\alpha \cdot T + \beta \cdot T_S) - k_{D_S} \cdot D_S $$ 
$$ \frac{dD_L}{dt} = \lambda_L \cdot (\alpha \cdot T + \beta \cdot T_S) - k_{D_L} \cdot D_L $$ 
$$ D_{S0} = D_S(t_0) $$ 
$$ D_{L0} = D_L(t_0) $$

230 where $t = 0$ denotes treatment initiation, and $T_{R0}$, $T_{S0}$, $D_{S0}$, and $D_{L0}$ are the initial observations of $T_R$, $T_S$, $D_S$, and $D_L$, respectively, at time $t_0$.

232 Statistical model
The parameters were identified using non-linear mixed effects modeling to account for inter-individual variability. Let $\mathcal{M}(t; \theta)$ be a structural model depending on time $t$ and a set $\theta$ of parameters. Then,

$$y_{ij} = \mathcal{M}(t_{ij}; \theta_i) + \epsilon_{ij}$$

where $y$ is observed data belonging to $\{T, D, D_S, D_L\}$, $y_{ij}$ is the j-th measurement of the i-th patient observed at time $t_{ij}$, $\theta_i \in \mathbb{R}^p$ is the set of p individual parameters of the patient $i$, and $\epsilon_{ij} \sim \mathcal{N}(\mu, \sigma_{ij})$ is the residual error model. Individual parameters $\theta_i = \{\theta_i^1, ..., \theta_i^p\}$ were assumed to follow a log-normal distribution, described by a fixed population parameter $\theta_{pop}$ and a random effect $\eta_i$ to account for inter-individual variability:

$$\ln(\theta_i) = \ln(\theta_{pop}) + \eta_i, \quad \eta_i \sim \mathcal{N}(0, \omega^2)$$

A constant error model was applied to $T$ ($\epsilon_{ij} \sim \mathcal{N}(1, a^2)$), whereas a proportional error model was applied to $D$, $D_S$, and $D_L$ ($\epsilon_{ij} \sim \mathcal{N}(0, (b \cdot \mathcal{M}(t_{ij}; \theta_i))^2)$).

The complete vector of population parameters $\psi$ was then extended to

$$\psi = (\theta_{pop}, \omega, a, b)$$

## Model calibration and selection

Population parameters $\psi$ of each structural model were identified by computing the maximum likelihood estimation using the Stochastic-Approximation Expectation Maximization algorithm$^{40}$ through Monolix software version 2023R1.

Tumor population parameters ($T_{S_0}$, $T_{R_0}$, $\alpha$, and $\beta$) were first identified from radiological report data, independently from cfDNA observations. They were secondly fixed in the complete joint cfDNA-TK model (Figure 1C).

Model was selected based on parameter identifiability (< 50% relative standard errors (RSEs) of population parameters $\psi$) and minimal corrected Bayesian Information Criterion (BICc) minimization$^{41}$:

$$\text{BICc} = -2\mathcal{L}_{\mathcal{M}}(\hat{\psi}) + \dim(\theta_R) \cdot \log(N) + \dim(\theta_F) \cdot \log\left(\sum_i n_i\right)$$

where $\mathcal{L}_{\mathcal{M}}(\hat{\psi})$ is the likelihood, $N$ is the number of patients, $\theta_R$ represents the number of random effects, and $\theta_F$ the number of fixed effects.
261 Individual empirical Bayes estimates were finally re-identified using only first 6 weeks
data from patients who had not progressed during this period, with fixed prior given by
262 the population parameters from the full kinetics model.
263

264 **Survival analysis**

265 To prevent immortal time bias due to model-based parameters calculated from post-
266 baseline data, we used landmark analysis⁴². We selected 6 weeks as a landmark time
267 point — corresponding to a minimum of two treatment cycles for most patients — and
268 defined PFS6 as the post-6 weeks PFS. PFS6 was computed as the time from 6 weeks
269 after treatment initiation to disease progression, death or last follow-up, whichever
270 occurred first. Progression date was defined as the date of disease progression based
271 on the iRECIST criteria.

272 All statistical calculations were performed using R version 4.4.0. Classical statistical
273 modeling was conducted for association with PFS6 using Cox proportional hazard
274 regression as implemented in the survival::coxph 3.8.3 function. Regressions were
275 adjusted for established prognostic clinical confounding variables: age, sex, ECOG
276 performance status, primary tumor origin and NLR. Multivariable Cox models including
277 either these clinical variables or these plus candidate additional biomarkers were
278 considered as predictive models and assessed by calculating the C-index using the
279 Hmisc::rcorr.cens 5.2.3 function. Survival was displayed using the Kaplan-Meier
280 estimated curves (survival::survfit 3.8.3). To optimally split patients into short and long
281 survivors based on continuous features, the function survminer::surv_cutpoint 0.5.0
282 was used, which maximizes the log-rank statistic. A minimum proportion of 20% of
283 patients per group was enforced.
284

# RESULTS

285

## Patient characteristics

286

287

288

289

290

291

292

293

The cohort included 112 patients suffering from advanced or metastatic carcinoma of multiple primary origin: UC (n = 8), HNSCC (n = 28), ccRCC (n = 13), or NSCLC (n = 63) (Table S1). Patients received either first- or second-line immunotherapy, most frequently pembrolizumab or nivolumab, with some in combination with chemotherapy or targeted therapy (Table S1). Pre-treatment SLD measurements were available for 95% of the cohort. The median PFS was 9.9 months [95% confidence interval (CI): 7.0-NA] months, with 52% of patients progressing for a median follow-up of 18.7 months (Table S1).

294

295

296

297

298

299

300

301

Tumor imaging was available with a median of 2 assessments per patient (range: 1-10), yielding a total of 295 tumor imaging, whereas cfDNA size profile quantification was performed with a median of 7 assessments per patient (range: 1-20), yielding a total of 907 cfDNA samples. Focusing on data available only at 6 weeks of treatment, tumor imaging was available once in 54 patients and twice in 39 patients. Quantification of the cfDNA size profile was available with a median of 3 assessments per patient (1 timepoint, n = 2; 2 timepoints: n = 25; 3 timepoints, n = 54; 4 timepoints, n = 12), with a total of 262 cfDNA samples.

302

## Tumor and cfDNA dynamic patterns

303

304

305

306

307

308

309

310

311

312

313

314

315

TK patterns varied widely among patients, mostly showing simple patterns of tumor growth or shrinkage (Figure 2A, Figure S1). Some patients also experienced an initial tumor shrinkage before subsequently progressing (n = 22). Two patients showed high tumor burden increase at treatment initiation (hyperprogression), and one patient demonstrated initial growth followed by treatment response (pseudoprogression). Total, short, and long cfDNA fragments also demonstrated high inter- and intra-patient variability over time. The concentration of all fragments ranged from 2.19 to 636.54 pg/µL, with a median of 13.47 pg/µL, (Figure 2B) samples >100 pg/µL. Twenty-three patients showed an initial spike in cfDNA levels at treatment initiation, followed by stabilization at lower levels, whereas other patients showed cfDNA spikes later during treatment (Figure S1). Most released fragments were approximately 166 bp and short fragments corresponded to 73% [range: 39-81%] of the total fragments, with a median concentration of 9.79 pg/µL [range: 1.38-512 pg/µL] (Figure 2C). Long fragments had
a smaller median concentration of 1.17 pg/µL [range: 0.19–44.63], with only 11 samples exceeding 10 pg/µL (Figure 2D). Moreover, long fragments did not follow the same temporal patterns as short fragments (Figure S1). Long fragments generally displayed concordant directional trends with short fragments, though with differing amplitudes; moreover, their dynamics occasionally diverged, in a way that increases in long fragments coincided with decreases in short fragments, or vice versa.

## 322 Mechanistic modeling of the joint tumor – size-dependent cfDNA kinetics

We developed a mechanistic model to jointly describe TK and size-structured short and long cfDNA fragment kinetics. The TK model and parameters were calibrated to the TK data.

The model was able to accurately describe non-trivial cfDNA kinetics (Figure 3A, Figure S2). Notably, it captured both steady trends and transient spikes at treatment initiation, observed in short and/or long fragment dynamics (e.g., patient 3, Figure 3A). It was also able to capture contrasting kinetic patterns between short and long fragments (e.g., a short fragment decrease with a long fragment increase in patient 2; or an initial short fragment decrease with an initial long fragment increase in patient 4).

Goodness-of-fit diagnostics indicated no model misspecification, supporting the adequacy of the model for jointly describing the observed tumor size dynamics to cfDNA concentration ones (Figure 3B-D). Individual conditionally weighted residuals of TK were centered around zero (Figure 3C-D), with an estimated residual error of 7.8 mm (Table 1A). However, the cfDNA ones displayed a small rightward skew, indicating a potential under-estimation of cfDNA data for all fragments, especially long ones (outlier proportion of 6.39% and 13% for short and long fragments, respectively, Figure 3B). Additionally, the population maximum likelihood estimation indicated a 45% and 57% discrepancy between the fits and the data for the short and long fragments, respectively (Table 1B).

All population parameter RSEs were below 19% and absolute correlations of the population parameter estimates did not exceed 0.31, indicating good parametric identifiability. A condition number of 3.46 further confirmed the model was not over-parameterized.

## 346 Quantitative insights into cfDNA biology
The mechanistic nature of the model allowed to derive quantitative parameters relative to the biology of cfDNA tumor release and elimination kinetics. Proliferative and dying cells were assumed to release cfDNA at the same shedding rate ($λ_S$ and $λ_L$ for short and long fragments, respectively), meaning that any change in the tumor cell population was converted into cfDNA release similarly. The shedding rate of small fragments was estimated to 0.46 pg.µL⁻¹.mm⁻¹, 7.4 times larger than for long fragments (Table 1B).

Important inter-patient variability was identified for $λ_S$ and $k_{D_S}$, with coefficients of variation of 240% and 250%, respectively (Table 1B), being both approximately 3 times higher than the variability of the baseline short fragment concentration. Similarly, $λ_L$ and $k_{D_L}$ showed coefficients of variation of 110% and 140%, respectively, being approximately 1.7 and 2.2-fold higher than the variability of the baseline long-fragment concentration.

We found moderate correlation of $λ_S$ with $λ_L$ ($r = 0.52, p < 0.0001$, Figure S4) and significant mean difference (Student's t-test, $p < 0.0001$). $k_{D_S}$ and $k_{D_L}$ showed moderate correlation ($r = 0.5, p < 0.001$) and significant difference (Student's t-test, $p < 0.001$).

**An early on-treatment model-derived parameter predicts progression-free survival beyond 6 weeks**

With population parameters held constant, individual parameters were re-estimated using data from the first six weeks of treatment for patients who had not progressed during this period (n = 93). Tumor kinetics parameters were associated with post-6 weeks PFS (PFS6, $\frac{\alpha}{\beta}$ C-Index = 0.64, HR: 1.6 [95% CI: 1.1-2.2], $p = 0.02$, Table 2, Figure 4AB), although with a wide confidence interval, reflecting the small number of patients with more than one measurement (n = 39). The ratio $\frac{\lambda_S}{k_{D_S}}$ was also significantly associated with shorter PFS6, even after adjusting for confounders, notably baseline NLR (HR: 1.6 [95% CI: 1.2 - 2.2], $p = 0.001$, Figure 4AB). Critically, this ratio provided predictive value beyond established prognostic variables. The multivariable Cox model, which achieved a C-index of 0.78 [95% CI: 0.73-0.89] when including only the clinical variables, significantly increased to 0.80 [95% CI: 0.74-0.90] after adding $\frac{\lambda_S}{k_{D_S}}$ (one-tailed paired t-test $p < 0.0001$).
378 Individual patient trajectories further illustrate the predictive potential of cfDNA-derived
379 model-based parameters (Figures 4C-S3). The optimal cutpoint for splitting PFS6
380 Kaplan-Meier curves (Figures 4B) was determined as $c_S = 0.36$ for $\frac{\lambda_S}{k_{DS}}$. This ratio was
381 above $c_S$ in patients 5 and 6, both of whom progressed within 3 months of treatment.
382 In contrast, patient 7, who progressed later at 8 months, and patient 8, who had not
383 progressed at 8 months, showed consistent values of $\frac{\lambda_S}{k_{DS}}$ below the identified cutpoints
384 (Figure 4C).

385
386 DISCUSSION

387 This study provides the first mechanistic modeling integrating tumor kinetics with size-
388 dependent cfDNA fragment dynamics under treatment (here, immunotherapy). The
389 framework was able to accurately describe both short and long fragment dynamics,
390 including non-trivial features such as transient cfDNA spikes at treatment initiation in
391 short and/or long fragments. The identifiability of model parameters enabled validation
392 of underlying biological assumptions, supporting a size-dependent release and
393 clearance of cfDNA fragments. Furthermore, a model-derived parameter – identifiable
394 for individual patients from early on-treatment data – was associated with PFS early
395 on-treatment, highlighting the model’s potential for therapeutic management.

396 The NLME mechanistic mathematical modeling approach provides multiple
397 advantages: 1) it allows determination of parameters integrating the entire kinetics
398 (instead of time-dependent snapshots); 2) these parameters have physiological
399 meaning; 3) the model fit allows smoothening of measurement errors, leading to
400 decreased variance, stronger robustness and, *in fine*, better predictions (as shown
401 previously²⁵) and 4) the NLME population approach makes it possible to derive
402 individual estimates from early on-treatment data through Bayesian estimation.

403 In this model, we distinguished between two size ranges of cfDNA fragments
404 (threshold at 580 bp) to capture the two main shedding processes: apoptosis and
405 necrosis. We considered the simplifying assumption that tumor cell death and
406 proliferation contributed equally to cfDNA release in the circulation. In a preliminary
407 version of the model, we had allowed for distinct shedding rates depending on whether
408 cells proliferated or died, but this resulting in parametric non-identifiability. Moreover,
409 the significant difference between k<sub>D<sub>S</sub></sub> and k<sub>D<sub>L</sub></sub>, along with their weak correlation,
410 indicates size-dependent clearance mechanisms, supporting previous results from
411 stochastic modeling of cfDNA fragmentation processes³⁰. Unexpectedly, we observed
412 a higher elimination rate for the long fragments (k<sub>D<sub>Spop</sub></sub> = 0.36 months⁻¹ versus k<sub>D<sub>Lpop</sub></sub> =
413 0.49 months⁻¹, Table 1). This contradicts previous reports that long fragments are less
414 efficiently cleared by kidneys⁴³, though it has only been proven for single-stranded
415 DNA. This discrepancy could be explained by the predominance of double-stranded
416 DNA, and nucleosome-bound fragments — mainly cleared by the liver⁴⁴ —, or by
417 additional fragmentation by nucleases within the circulation, causing the concentration
418 of long fragments to decrease while the one of short fragments to increase. In early
419 versions of the model, we considered incorporating a fragmentation process in which
420 long fragments break down into short ones. This approach was ultimately discarded,
421 because no effective fragmentation term was parametrically identifiable on the time
422 scale of the available data (every 2-3 weeks), as longer fragments are continuously
423 degraded in the circulation by nuclease activity. Collecting cfDNA data at a finer time
424 scale could enable integration and better characterization of such a fragmentation
425 process.

426 Short fragments release during the first six weeks of treatment was significantly
427 associated with post-6 weeks PFS. These fragments comprise the first three peaks in
428 cfDNA size profiles and correspond to mono-, di-, and tri-nucleosomal fragments and
429 are primarily released through apoptosis³. This form of programmed cell death process
430 involves cfDNA fragmentation through nuclease activation within cell nucleus, which
431 cleaves nucleotides around histones and leads to these characteristic patterns. This
432 association likely reflects how treatment-sensitive tumor cells respond to
433 immunotherapy, indicating that short fragment dynamics at treatment initiation may
434 serve as a systemic surrogate marker of response by quantifying the proportion of cells
435 undergoing apoptosis. Crucially, this information is readily obtainable from plasma
436 samples, offering a cost-effective, accessible, and frequently measurable approach to
437 monitor emerging treatment resistance. Notably, this parameter could only be derived
438 through mechanistic modeling and no predictive kinetic parameter was found from the
439 raw data or through simpler, empirical, modeling approaches. In clinical practice, this
440 tool could inform adaptive therapeutic strategies, such as accelerating imaging-based
441 response assessment or, in conjunction with other factors like treatment intolerance,
442 redirecting patients toward alternative therapies.

443 Despite these strengths, our work also presents some limitations. This study only
444 considered general tumor growth and shrinkage, by calibrating the model on the SLD
445 of all lesions, without considering individual lesion kinetics. Future work could go
446 further on the modeling of individual lesions following recent advanced by Kerioui et
447 al.⁴⁵,⁴⁶ and, possibly, organ-specific shedding rates. Conversely, the current model still
448 relies on (limited) tumor size imaging assessments. Further research will explore the
449 predictive value of a model relying on cfDNA data only. Moreover, cfDNA release was
450 modeled as tumor-driven only, whereas white blood cells have been shown to release
from 55% of the cfDNA in healthy individuals (32% granulocytes, 11% monocytes and 12% lymphocytes)⁴⁷ to 76% in cancer patients⁴⁸. This omission may explain the overestimation of long-fragment spikes observed later during treatment in some patients (e.g., patient 7 at ~3 months, Figure 3). Integrating immune cell kinetics (e.g., neutrophils, lymphocytes) would surely improve long cfDNA fragment estimation. Methodologically, early on-treatment parameter estimation was performed on the same cohort used for population calibration, restricted to patients who did not progress during the considered period. This highlights the need for independent validation in future cohorts. Moreover, PFS has been considered the primary endpoint, but additional outcomes could be explored in the future, such as OS. More advanced statistical approaches, such as joint modeling⁴⁹ could be explored to improve accuracy and account for informative dropout (early study exit because of early progression or death). Finally, only univariable predictive models have been explored here. Multivariable analysis, integrating pre-treatment and early on-treatment model-based parameters with clinical and biological covariates, could be further explored, using, e.g., machine learning frameworks⁷. These could improve treatment response prediction and go deeper in the use of cfDNA as biomarker of treatment monitoring.
Altogether, this study demonstrates that mechanistic modeling of size-dependent cfDNA fragments provides biological insights and early predictive markers of ICI response, highlighting the causal role of tumor dynamics on cfDNA concentrations by size. Beyond its immediate clinical applications, this approach opens the way for digital twin models in oncology, enabling in silico testing of hypotheses about cfDNA release and clearance, e.g., simulating the effect of therapeutic interventions targeting the tumor on cfDNA kinetics, and potentially guiding clinical trial design. With further validation and integration of immune system dynamics, mechanistic cfDNA models could become powerful tools for precision immuno-oncology.

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

480
481 (A) cfDNA size profile of one patient prior to treatment. Both x- and y-axes are displayed on logarithmic scales. The concentration (area under the curve) of short 482 fragments (75-580 base pairs (bp)) is shown in green, while the concentration of long 483 fragments (580-1650 bp) is shown in blue. (B) Joint cfDNA-TK model schematics. 484
485 Tumor cells proliferate at a rate of $\alpha$. Treatment-sensitive tumor cells die at a rate of $\beta$. 486 cfDNA fragments are categorized as either short or long. Both fragment types are 487 released by proliferating and dying cells, at size-specific rates ($\lambda_S$ and $\lambda_L$, respectively). 488
489 Fragments are proportionally cleared at size-dependent rates ($k_{DS}$ and $k_{DL}$). (C) Sequential parameter estimation framework.
<figure><img src="image_3.png" /><figcaption>**Figure 2: Spaghetti plots of tumor and cfDNA kinetics**</figcaption></figure>

Individual patient tumor size and cfDNA concentration kinetics over time. (A) Tumor size kinetics for all 112 patients of the cohort, represented by the sum of largest diameters (SLD) of all lesions, in millimeters (mm). (B-D) Concentration (area under the size profile) of all fragments (75-1650 base pairs (bp), B), short fragments (75-580 bp, C), and long fragments (580-1650 bp, D). All panels use a logarithmic scale on the y-axis.
<figure><img src="image_4.png" /><figcaption>**Figure 3: The joint size-dependent cfDNA-tumor kinetics model accurately fits size-dependent cfDNA kinetics**</figcaption></figure>
(A) Individual model fits for 4 patients using the joint size-dependent cfDNA-TK model. Top panels: variable $T'$, representing tumor size kinetics, in red. Data points are the observed sum of largest diameters (SLD) of all lesions, in millimeters (mm). Middle panels: variable $D_S$, representing short fragment kinetics, in green. Data points are the observed cfDNA concentration of short fragments (75-580 base pairs (bp)). Bottom panels: variable $D_L$, representing long fragment kinetics, in blue. Data points are the observed cfDNA concentration of short fragments (580-1650 bp). (B-D) Goodness-of-fit diagnostic plot for tumor (red), short fragment (green), and long fragment (blue) kinetics. (B) Individual predictions vs. observations. Outliers are defined as predictions outside of the 90% prediction interval. (C) Individual predictions vs. individual weighted residuals (IWRES). $\epsilon$-shrinkage represents the information on residual variability and is defined as $(1 - \sigma_{IWRES}) \times 100$, where $\sigma_{IWRES}$ is the standard deviation of the IWRES$^50$. Informative data lead to an $\epsilon$-shrinkage close to 0%, whereas non-informative data will result in $\epsilon$-shrinkage values approaching 100%. (D) Empirical IWRES distribution. Dotted lines represent the theoretical probability density function.
<figure><img src="image_5.png" /><figcaption>**Figure 4: Short fragment-based parameters estimated early during treatment are associated with longer PFS**</figcaption></figure>

(A) Cox proportional hazards regression results for the best model-derived parameters estimated from early on-treatment kinetics (within 6 weeks of therapy), displayed through forest plot of the hazard ratios and their 95% confidence interval (CI). Corresponding C-index are displayed on the right of the forest plots. (B) Kaplan-Meier curves of progression-free survival stratified by model parameters. The log-rank test p-value and the optimal threshold used to separate patients with long vs. short progression-free survival (see Methods) are shown. (C) Early on-treatment individual
model fits for four patients using the joint size-dependent cfDNA-TK model, estimated
using data available at 6 weeks of treatment. Patients 5 and 6 progressed before 3
months, while patient 7 progressed at 8 months and patient 8 did not progress. Top
panels: variable *T*, representing tumor size kinetics, in red. Data points are the
observed sum of largest diameters (SLD) of all lesions, in millimeters (mm). Middle
panels: variable *D*<sub>s</sub>, representing short fragment kinetics, in green. Data points are the
observed cfDNA concentration of short fragments (75-580 base pairs (bp)). Bottom
panels: variable *D*<sub>L</sub>, representing long fragment kinetics, in blue. Data points are the
observed cfDNA concentration of short fragments (580-1650 bp).
TABLES

<table><caption>Table 1: Population parameters of the joint size-dependent cfDNA-TK model</caption><thead><tr><th colspan="5">TK model</th><th colspan="5">Joint model</th></tr><tr><th></th><th>UNIT</th><th>VALUE</th><th>C.V.(%)</th><th>STOCHASTIC<br>APPROXIMATION</th><th></th><th>UNIT</th><th>VALUE</th><th>C.V.(%)</th><th>STOCHASTIC<br>APPROXIMATION</th></tr><tr><th></th><th></th><th></th><th></th><th>R.S.E.</th><th></th><th></th><th></th><th></th><th>R.S.E.</th></tr><tr><th></th><th></th><th></th><th></th><th>(%)</th><th></th><th></th><th></th><th></th><th>(%)</th></tr></thead><tbody><tr><th colspan="5">FIXED EFFECTS</th><th colspan="5">FIXED EFFECTS</th></tr><tr><td>α<sub>pop</sub></td><td>month<sup>-1</sup></td><td>0.12</td><td></td><td>0.019</td><td>16</td><td>D<sub>S<sub>0</sub>pop</sub></td><td>pg·μL<sup>-1</sup></td><td>11</td><td>0.99</td><td>9.2</td></tr><tr><td>β<sub>pop</sub></td><td>month<sup>-1</sup></td><td>0.38</td><td></td><td>0.043</td><td>11</td><td>D<sub>L<sub>0</sub>pop</sub></td><td>pg·μL<sup>-1</sup></td><td>1.3</td><td>0.11</td><td>8.3</td></tr><tr><td>T<sub>R<sub>0</sub>pop</sub></td><td>mm</td><td>7.8</td><td></td><td>1.5</td><td>19</td><td>λ<sub>S<sub>pop</sub></sub></td><td>pg·μL<sup>-1</sup>·mm<sup>-1</sup></td><td>0.46</td><td>0.074</td><td>16</td></tr><tr><td>T<sub>S<sub>0</sub>pop</sub></td><td>mm</td><td>48</td><td></td><td>4.5</td><td>9.4</td><td>λ<sub>L<sub>pop</sub></sub></td><td>pg·μL<sup>-1</sup>·mm<sup>-1</sup></td><td>0.062</td><td>0.0075</td><td>12</td></tr><tr><th colspan="5">STANDARD DEVIATION OF THE RANDOM EFFECTS</th><th colspan="5">STANDARD DEVIATION OF THE RANDOM EFFECTS</th></tr><tr><td>ω<sub>α</sub></td><td></td><td>0.87</td><td>110</td><td>0.12</td><td>13</td><td>ω<sub>D<sub>S<sub>0</sub></sub></sub></td><td></td><td>0.68</td><td>77</td><td>0.078</td><td>11</td></tr><tr><td>ω<sub>β</sub></td><td></td><td>0.79</td><td>93</td><td>0.11</td><td>14</td><td>ω<sub>D<sub>L<sub>0</sub></sub></sub></td><td></td><td>0.59</td><td>65</td><td>0.069</td><td>12</td></tr><tr><td>ω<sub>T<sub>R<sub>0</sub></sub></sub></td><td></td><td>0.96</td><td>120</td><td>0.16</td><td>17</td><td>ω<sub>λ<sub>S</sub></sub></td><td></td><td>1.4</td><td>240</td><td>0.12</td><td>9.1</td></tr><tr><td>ω<sub>T<sub>S<sub>0</sub></sub></sub></td><td></td><td>0.88</td><td>110</td><td>0.069</td><td>7.9</td><td>ω<sub>λ<sub>L</sub></sub></td><td></td><td>0.92</td><td>110</td><td>0.075</td><td>8.2</td></tr><tr><th colspan="5">ERROR MODEL PARAMETERS</th><th colspan="5">ERROR MODEL PARAMETERS</th></tr><tr><td>a</td><td>mm</td><td>7.8</td><td></td><td>0.59</td><td>7.6</td><td>ω<sub>k<sub>D<sub>S</sub></sub></sub></td><td></td><td>1.4</td><td>250</td><td>0.13</td><td>9.4</td></tr><tr><td></td><td></td><td></td><td></td><td></td><td></td><td>ω<sub>k<sub>D<sub>L</sub></sub></sub></td><td></td><td>1.1</td><td>140</td><td>0.1</td><td>9.7</td></tr><tr><th colspan="5">ERROR MODEL PARAMETERS</th><th colspan="5">ERROR MODEL PARAMETERS</th></tr><tr><td>b<sub>SHORT</sub></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td>0.45</td><td></td><td>0.014</td><td>3.1</td></tr><tr><td>a<sub>TK</sub></td><td>mm</td><td></td><td></td><td></td><td></td><td></td><td></td><td>7.8</td><td></td><td></td><td></td></tr><tr><td>b<sub>LONG</sub></td><td></td><td></td><td></td><td></td><td></td><td></td><td></td><td>0.57</td><td></td><td>0.019</td><td>3.3</td></tr></tbody></table>

(C.V.: coefficient of variation, defined as $100 \cdot \sqrt{\exp(\omega_{\theta}^2) - 1}$. S.E.: standard error.

R.S.E.: relative standard error. TK: tumor kinetic.)
<table><caption>Table 2: Cox proportional hazards regression results</caption><thead><tr><th rowspan="3">VARIABLE</th><th rowspan="3">LEVEL</th><th rowspan="3">N</th><th colspan="3">EARLY KINETICS (6 WEEKS)</th></tr><tr><th colspan="3">PROGRESSION-FREE SURVIVAL</th></tr><tr><th>C INDEX</th><th>HR</th><th>SIGNIF</th></tr></thead><tbody><tr><td colspan="6">MODEL-BASED PARAMETERS</td></tr><tr><td>α/β</td><td></td><td>93</td><td>0.64</td><td>1.6 (1.1 – 2.2)</td><td>**</td></tr><tr><td>λ<sub>S</sub>/k<sub>D<sub>S</sub></sub></td><td></td><td>93</td><td>0.59</td><td>1.6 (1.2 – 2.2)</td><td>**</td></tr><tr><td>λ<sub>S</sub></td><td></td><td>93</td><td>0.57</td><td>1.6 (1.2 – 2.1)</td><td>***</td></tr><tr><td colspan="6">CLINICAL AND BIOLOGICAL VARIABLES</td></tr><tr><td>Pathology</td><td>NSCLC (reference)</td><td>93</td><td></td><td></td><td></td></tr><tr><td></td><td>HNSCC</td><td></td><td>0.66</td><td>6.6 (2.8 – 15)</td><td>****</td></tr><tr><td></td><td>ccRCC</td><td></td><td>0.66</td><td>2 (0.74 – 5.6)</td><td>n.s.</td></tr><tr><td></td><td>UC</td><td></td><td>0.66</td><td>3.6 (0.88 – 15)</td><td>n.s.</td></tr><tr><td>NLR</td><td></td><td>84</td><td>0.63</td><td>1.2 (0.82 – 1.9)</td><td>n.s.</td></tr><tr><td>Sex</td><td>female (reference)</td><td>93</td><td></td><td></td><td></td></tr><tr><td></td><td>male</td><td></td><td>0.61</td><td>4.6 (1.6 – 14)</td><td>**</td></tr><tr><td>ECOG</td><td>0 (reference)</td><td>93</td><td></td><td></td><td></td></tr><tr><td></td><td>2+</td><td></td><td>0.56</td><td>1.5 (0.58 – 3.8)</td><td>n.s.</td></tr><tr><td></td><td>1</td><td></td><td>0.56</td><td>0.6 (0.27 – 1.4)</td><td>n.s.</td></tr><tr><td>Age</td><td></td><td>93</td><td>0.50</td><td>0.88 (0.61 – 1.3)</td><td>n.s.</td></tr></tbody></table>

Results are displayed for the model-based parameters derived from early on-treatment kinetics (within 6 weeks of therapy) and the clinical pre-treatment variables.

(ccRCC: clear cell renal cell carcinoma; HNSCC: head and neck squamous cell carcinoma; HR: hazard ratio; MV: multivariable; NSCLC: non-small cell lung cancer; OR: odds ratio; UC: urothelial carcinoma; UV: univariable)

(Significance: ****: p-value < 0.0001; ***: p-value < 0.001; **: p-value < 0.01; *: p-value < 0.05; n.s.: non-significant, > 0.05)
SUPPLEMENTARY

<caption>**Figure S1: Exploratory longitudinal data analysis**</caption>

See PDF

<caption>**Figure S2: Individual fits – estimation from full kinetics**</caption>

See PDF

<caption>**Figure S3: Individual fits – estimation from early kinetics**</caption>

See PDF

<figure><img src="image_6.png" /><figcaption>**Table s1: Patient characteristics**</figcaption></figure>
<table><thead><tr><th rowspan="2">Variable</th><th colspan="5">Pathology</th><th rowspan="2">p-value<sup>2</sup></th></tr><tr><th>Overall, N = 112<sup>1</sup></th><th>UC, N = 8<sup>1</sup></th><th>HNSCC, N = 28<sup>1</sup></th><th>ccRCC, N = 13<sup>1</sup></th><th>NSCLC, N = 63<sup>1</sup></th></tr></thead><tbody><tr><td>Sex</td><td></td><td></td><td></td><td></td><td></td><td>0,9</td></tr><tr><td>Female</td><td>31 (28%)</td><td>1 (13%)</td><td>8 (29%)</td><td>4 (31%)</td><td>18 (29%)</td><td></td></tr><tr><td>Male</td><td>81 (72%)</td><td>7 (88%)</td><td>20 (71%)</td><td>9 (69%)</td><td>45 (71%)</td><td></td></tr><tr><td>Age (years)</td><td>64 (11)</td><td>71 (9)</td><td>66 (9)</td><td>58 (14)</td><td>63 (11)</td><td>0,10</td></tr><tr><td>Treatment<sup>3</sup></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Atezolizumab</td><td>1 (0,9%)</td><td>0 (0%)</td><td>0 (0%)</td><td>0 (0%)</td><td>1 (1,6%)</td><td></td></tr><tr><td>Nivolumab</td><td>26 (23%)</td><td>0 (0%)</td><td>20 (71%)</td><td>0 (0%)</td><td>6 (9,5%)</td><td></td></tr><tr><td>Nivolumab + Ipilimumab</td><td>2 (1,8%)</td><td>0 (0%)</td><td>0 (0%)</td><td>2 (15%)</td><td>0 (0%)</td><td></td></tr><tr><td>Nivolumab + targeted therapy</td><td>4 (3,6%)</td><td>0 (0%)</td><td>0 (0%)</td><td>4 (31%)</td><td>0 (0%)</td><td></td></tr><tr><td>Pembrolizumab</td><td>27 (24%)</td><td>7 (88%)</td><td>2 (7,1%)</td><td>2 (15%)</td><td>16 (25%)</td><td></td></tr><tr><td>Pembrolizumab + chemotherapy</td><td>45 (40%)</td><td>0 (0%)</td><td>5 (18%)</td><td>0 (0%)</td><td>40 (63%)</td><td></td></tr><tr><td>Pembrolizumab + targeted therapy</td><td>7 (6,3%)</td><td>1 (13%)</td><td>1 (3,6%)</td><td>5 (38%)</td><td>0 (0%)</td><td></td></tr><tr><td>Treatment line</td><td></td><td></td><td></td><td></td><td></td><td>&lt;0,001</td></tr><tr><td>1</td><td>91 (81%)</td><td>3 (38%)</td><td>19 (68%)</td><td>13 (100%)</td><td>56 (89%)</td><td></td></tr><tr><td>2</td><td>21 (19%)</td><td>5 (63%)</td><td>9 (32%)</td><td>0 (0%)</td><td>7 (11%)</td><td></td></tr><tr><td>CPS (%)</td><td></td><td></td><td></td><td></td><td></td><td>0,6</td></tr><tr><td>[1, 20[</td><td>7 (39%)</td><td>0 (0%)</td><td>7 (41%)</td><td>0 (-%)</td><td>0 (-%)</td><td></td></tr><tr><td>&lt; 1</td><td>5 (28%)</td><td>0 (0%)</td><td>5 (29%)</td><td>0 (-%)</td><td>0 (-%)</td><td></td></tr><tr><td>&gt;= 20</td><td>6 (33%)</td><td>1 (100%)</td><td>5 (29%)</td><td>0 (-%)</td><td>0 (-%)</td><td></td></tr><tr><td>Missing</td><td>94</td><td>7</td><td>11</td><td>13</td><td>63</td><td></td></tr><tr><td>TPS (%)</td><td></td><td></td><td></td><td></td><td></td><td>&gt;0,9</td></tr><tr><td>[1, 50[</td><td>20 (35%)</td><td>0 (-%)</td><td>0 (-%)</td><td>0 (-%)</td><td>20 (35%)</td><td></td></tr><tr><td>&lt; 1</td><td>14 (25%)</td><td>0 (-%)</td><td>0 (-%)</td><td>0 (-%)</td><td>14 (25%)</td><td></td></tr><tr><td>&gt;= 50</td><td>23 (40%)</td><td>0 (-%)</td><td>0 (-%)</td><td>0 (-%)</td><td>23 (40%)</td><td></td></tr><tr><td>Missing</td><td>55</td><td>8</td><td>28</td><td>13</td><td>6</td><td></td></tr><tr><td>NLR</td><td>4,5 (3,7)</td><td>4,2 (2,2)</td><td>7,3 (5,3)</td><td>2,7 (1,2)</td><td>3,6 (2,3)</td><td>0,004</td></tr><tr><td>Missing</td><td>10</td><td>0</td><td>2</td><td>1</td><td>7</td><td></td></tr><tr><td>SLD (mm)</td><td>67 (50)</td><td>49 (42)</td><td>39 (23)</td><td>99 (59)</td><td>74 (52)</td><td>0,001</td></tr><tr><td>Missing</td><td>3</td><td>0</td><td>1</td><td>0</td><td>2</td><td></td></tr><tr><td>Mutations<sup>4</sup></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>KRAS</td><td>24 (63%)</td><td>0 (0%)</td><td>0 (0%)</td><td>0 (-%)</td><td>24 (67%)</td><td></td></tr><tr><td>TP53</td><td>24 (71%)</td><td>1 (100%)</td><td>1 (100%)</td><td>0 (-%)</td><td>22 (69%)</td><td></td></tr><tr><td>Early progression<sup>5</sup></td><td></td><td></td><td></td><td></td><td></td><td>&lt;0,001</td></tr><tr><td>0</td><td>76 (68%)</td><td>3 (38%)</td><td>12 (43%)</td><td>10 (77%)</td><td>51 (81%)</td><td></td></tr><tr><td>1</td><td>36 (32%)</td><td>5 (63%)</td><td>16 (57%)</td><td>3 (23%)</td><td>12 (19%)</td><td></td></tr><tr><td>PFS (median)</td><td>9,9 (7,0, -)</td><td>2,9 (2,4, 6,8)</td><td>2,6 (1,9, -)</td><td>9,6 (5,8, -)</td><td>- (13, -)</td><td>&lt;0,001</td></tr></tbody></table>

(ccRCC: clear cell renal cell carcinoma; CPS: combined positive score; HNSCC: head and neck squamous cell carcinoma; NLR: neutrophil-to-lymphocyte ratio; NSCLC: non-small cell lung cancer; PFS: progression-free survival; SLD: sum of the longest diameter of the lesions; TPS: tumor proportion score; UC: urothelial carcinoma;)

<footnote n="1">¹ n (%); Mean (SD)</footnote>

<footnote n="2">² Fisher's exact test; Kruskal-Wallis rank sum test; NA</footnote>

<footnote n="3">³ Immunotherapy combination and potential associated therapy</footnote>

<footnote n="4">⁴ Other hidden mutations are: EGFR (3), BRAF (4), CTNNB1 (1), ERBB4 (1), FGFR 3 (1), P13 KINASE (1), PDGFRA (1), PIK3CA (2), PTEN (4), RET (1), STK11 (1), TTF1 (1)</footnote>

<footnote n="5">⁵ 1, progression before 3 months of immunotherapy | 0, otherwise</footnote>
# ACKNOWLEDGEMENTS

This work received support from the French government under the France 2030 investment plan, as part of the Initiative d'Excellence d'Aix-Marseille Université - A*MIDEX (AMX-19-IET-001 & AMX-21-IET-017).

# Declaration of Interest statement

FF, FG and AB are all members of Adelis Technologies, Grabels, France, which sells the BIABooster platform. FF is member of IDS-Solutions Oncology (Marseille, France), which is shareholder of Adelis Technologies.

# REFERENCES

1. Ivanov M, Baranova A, Butler T, Spellman P, Mileyko V. Non-random fragmentation patterns in circulating cell-free DNA reflect epigenetic regulation. *BMC Genom*. 2015;16(Suppl 13):S1. doi:10.1186/1471-2164-16-S13-S1

2. Qi T, Pan M, Shi H, Wang L, Bai Y, Ge Q. Cell-Free DNA Fragmentomics: The Novel Promising Biomarker. *Int J Mol Sci*. 2023;24(2):1503. doi:10.3390/ijms24021503

3. Heitzer E, Auinger L, Speicher MR. Cell-Free DNA and Apoptosis: How Dead Cells Inform About the Living. *Trends in Molecular Medicine*. 2020;26(5):519-528. doi:10.1016/j.molmed.2020.01.012

4. Sharma P, Hu-Lieskovan S, Wargo JA, Ribas A. Primary, Adaptive, and Acquired Resistance to Cancer Immunotherapy. *Cell*. 2017;168(4):707-723. doi:10.1016/j.cell.2017.01.017

5. Seymour L, Bogaerts J, Perrone A, et al. iRECIST: guidelines for response criteria for use in trials testing immunotherapeutics. *Lancet Oncol*. 2017;18(3):e143-e152. doi:10.1016/S1470-2045(17)30074-8

6. Lu S, Stein JE, Rimm DL, et al. Comparison of Biomarker Modalities for Predicting Response to PD-1/PD-L1 Checkpoint Blockade: A Systematic Review and Meta-analysis. *JAMA Oncology*. 2019;5(8):1195-1204. doi:10.1001/jamaoncol.2019.1549

7. Nguyen Phuong L, Salas S, Benzekry S. Computational Modeling for Circulating Cell-Free DNA in Clinical Oncology. *JCO Clin Cancer Inform*. 2025;(9):e2400224. doi:10.1200/CCI-24-00224
8. Si H, Kuziora M, Quinn KJ, et al. A Blood-based Assay for Assessment of Tumor Mutational Burden in First-line Metastatic NSCLC Treatment: Results from the MYSTIC Study. *Clinical Cancer Research*. 2021;27(6):1631-1640. doi:10.1158/1078-0432.CCR-20-3771

9. Gandara DR, Paul SM, Kowanetz M, et al. Blood-based tumor mutational burden as a predictor of clinical benefit in non-small-cell lung cancer patients treated with atezolizumab. *Nat Med*. 2018;24(9):1441-1448. doi:10.1038/s41591-018-0134-3

10. Jensen TJ, Goodman AM, Kato S, et al. Genome-Wide Sequencing of Cell-Free DNA Identifies Copy-Number Alterations That Can Be Used for Monitoring Response to Immunotherapy in Cancer Patients. *Molecular Cancer Therapeutics*. 2019;18(2):448-458. doi:10.1158/1535-7163.MCT-18-0535

11. Vargas-Accarino E, Higuera M, Bermúdez-Ramos M, et al. Harnessing Plasma Biomarkers to Predict Immunotherapy Outcomes in Hepatocellular Carcinoma: The Role of cfDNA, ctDNA, and Cytokines. *International Journal of Molecular Sciences*. 2025;26(6):2794. doi:10.3390/ijms26062794

12. Lapin M, Oltedal S, Tjensvoll K, et al. Fragment size and level of cell-free DNA provide prognostic information in patients with advanced pancreatic cancer. *J Transl Med*. 2018;16:300. doi:10.1186/s12967-018-1677-2

13. van 't Erve I, Alipanahi B, Lumbard K, et al. Cancer treatment monitoring using cell-free DNA fragmentomes. *Nat Commun*. 2024;15(1):8801. doi:10.1038/s41467-024-53017-7

14. Janke F, Gasser M, Angeles AK, et al. Low-coverage whole genome sequencing of cell-free DNA to predict and track immunotherapy response in advanced non-small cell lung cancer. *Journal of Experimental & Clinical Cancer Research*. 2025;44(1):87. doi:10.1186/s13046-025-03348-0

15. Bratman SV, Yang SYC, Iafolla MAJ, et al. Personalized circulating tumor DNA analysis as a predictive biomarker in solid tumor patients treated with pembrolizumab. *Nat Cancer*. 2020;1(9):873-881. doi:10.1038/s43018-020-0096-5

16. Moding EJ, Liu Y, Nabet BY, et al. Circulating Tumor DNA Dynamics Predict Benefit from Consolidation Immunotherapy in Locally Advanced Non-Small Cell Lung Cancer. *Nat Cancer*. 2020;1(2):176-183. doi:10.1038/s43018-019-0011-0

17. Lipson EJ, Velculescu VE, Pritchard TS, et al. Circulating tumor DNA analysis as a real-time method for monitoring tumor burden in melanoma patients undergoing treatment with
immune checkpoint blockade. *J Immunother Cancer*. 2014;2(1):42. doi:10.1186/s40425-014-0042-0

18. Cabel L, Riva F, Servois V, et al. Circulating tumor DNA changes for early monitoring of anti-PD1 immunotherapy: a proof-of-concept study. *Ann Oncol*. 2017;28(8):1996-2001. doi:10.1093/annonc/mdx212

19. Assaf ZJF, Zou W, Fine AD, et al. A longitudinal circulating tumor DNA-based model associated with survival in metastatic non-small-cell lung cancer. *Nat Med*. Published online March 16, 2023:1-10. doi:10.1038/s41591-023-02226-6

20. Ribba B, Holford N, Magni P, et al. A Review of Mixed-Effects Models of Tumor Growth and Effects of Anticancer Drug Treatment Used in Population Analysis. *CPT: Pharmacometrics & Systems Pharmacology*. 2014;3(5):113. doi:10.1038/psp.2014.12

21. Sheiner LB. The Population Approach to Pharmacokinetic Data Analysis: Rationale and Standard Data Analysis Methods. *Drug Metabolism Reviews*. 1984;15(1-2):153-171. doi:10.3109/03602538409015063

22. Almufti R, Wilbaux M, Oza A, et al. A critical review of the analytical approaches for circulating tumor biomarker kinetics during treatment. *Annals of Oncology*. 2014;25(1):41-56. doi:10.1093/annonc/mdt382

23. Hansson E, Amantea M, Westwood P, et al. PKPD Modeling of VEGF, sVEGFR-2, sVEGFR-3, and sKIT as Predictors of Tumor Dynamics and Overall Survival Following Sunitinib Treatment in GIST. *CPT: Pharmacometrics & Systems Pharmacology*. 2013;2(11):84. doi:10.1038/psp.2013.61

24. Netterberg I, Karlsson MO, Terstappen LWMM, Koopman M, Punt CJA, Friberg LE. Comparing Circulating Tumor Cell Counts with Dynamic Tumor Size Changes as Predictor of Overall Survival: A Quantitative Modeling Framework. *Clin Cancer Res*. 2020;26(18):4892-4900. doi:10.1158/1078-0432.CCR-19-2570

25. Benzekry S, Karlsen M, Bigarré C, et al. Predicting Survival in Patients with Advanced NSCLC Treated with Atezolizumab Using Pre- and on-Treatment Prognostic Biomarkers. *Clinical Pharmacology & Therapeutics*. 2024;116(4):1110-1120. doi:10.1002/cpt.3371

26. Avanzini S, Kurtz DM, Chabon JJ, et al. A mathematical model of ctDNA shedding predicts tumor detection size. *Sci Adv*. 2020;6(50):eabc4308. doi:10.1126/sciadv.abc4308

27. Janssen JM, Verheijen RB, van Duijl TT, et al. Longitudinal nonlinear mixed effects modeling of EGFR mutations in ctDNA as predictor of disease progression in treatment of
EGFR-mutant non-small cell lung cancer. Clin Transl Sci. 2022;15(8):1916-1925.
doi:10.1111/cts.13300

28. Li A, Lou E, Leder K, Foo J. Early Circulating Tumor DNA Kinetics as a Dynamic Biomarker of Cancer Treatment Response. JCO Clin Cancer Inform. 2025;9:e2400160.
doi:10.1200/CCI-24-00160

29. Ribba B, Roller A, Helms HJ, Stern M, Bleul C. Circulating tumor DNA: Opportunities and challenges for pharmacometric approaches. Front Pharmacol. 2023;13. Accessed October 26, 2023. https://www.frontiersin.org/articles/10.3389/fphar.2022.1058220

30. Tsui THL, Xie PF, Chulián S, Pérez-García VM. Breaking Down Cell-Free DNA Fragmentation: A Markov Model Approach. bioRxiv. Preprint posted online July 7, 2023;2023.07.06.547953. doi:10.1101/2023.07.06.547953

31. Jahr S, Hentze H, Englisch S, et al. DNA Fragments in the Blood Plasma of Cancer Patients: Quantitations and Evidence for Their Origin from Apoptotic and Necrotic Cells1. Cancer Res. 2001;61(4):1659-1665.

32. Nguyen Phuong L, Fina F, Greillier L, et al. The SChISM study: Cell-free DNA size profiles as predictors of progression in advanced carcinoma treated with immune-checkpoint inhibitors. Published online September 3, 2025. Accessed September 16, 2025. https://inria.hal.science/hal-05238567

33. Boutonnet A, Pradines A, Mano M, et al. Size and Concentration of Cell-Free DNA Measured Directly from Blood Plasma, without Prior DNA Extraction. Anal Chem. 2023;95(24):9263-9270. doi:10.1021/acs.analchem.3c00998

34. Bagley SJ, Kothari S, Aggarwal C, et al. Pretreatment neutrophil-to-lymphocyte ratio as a marker of outcomes in nivolumab-treated patients with advanced non-small-cell lung cancer. Lung Cancer. 2017;106:1-7. doi:10.1016/j.lungcan.2017.01.013

35. Mezquita L, Auclin E, Ferrara R, et al. Association of the Lung Immune Prognostic Index With Immune Checkpoint Inhibitor Outcomes in Patients With Advanced Non-Small Cell Lung Cancer. JAMA Oncology. 2018;4(3):351-357. doi:10.1001/jamaoncol.2017.4771

36. Khier S, Gahan PB. Hepatic Clearance of Cell-Free DNA: Possible Impact on Early Metastasis Diagnosis. Mol Diagn Ther. 2021;25(6):677-682. doi:10.1007/s40291-021-00554-2

37. Kustanovich A, Schwartz R, Peretz T, Grinshpun A. Life and death of circulating cell-free DNA. Cancer Biol Ther. 2019;20(8):1057-1067. doi:10.1080/15384047.2019.1598759
38. Grabuschnig S, Bronkhorst AJ, Holdenrieder S, et al. Putative Origins of Cell-Free DNA in Humans: A Review of Active and Passive Nucleic Acid Release Mechanisms. *Int J Mol Sci.* 2020;21(21):8062. doi:10.3390/ijms21218062

39. Bronkhorst AJ, Wentzel JF, Ungerer V, et al. Sequence analysis of cell-free DNA derived from cultured human bone osteosarcoma (143B) cells. *Tumour Biol.* 2018;40(9):1010428318801190. doi:10.1177/1010428318801190

40. Kuhn E, Lavielle M. Maximum likelihood estimation in nonlinear mixed effects models. *Computational Statistics & Data Analysis.* 2005;49(4):1020-1038. doi:10.1016/j.csda.2004.07.002

41. Delattre M, Lavielle M, Poursat MA. A note on BIC in mixed-effects models. *Electronic Journal of Statistics*. 2014;8:456-475. doi:10.1214/14-EJS890

42. Khandelwal A, Grisic AM, French J, Venkatakrishnan K. Pharmacometrics Golems: Exposure-Response Models in Oncology. *Clinical Pharmacology & Therapeutics.* 2022;112(5):941-945. doi:10.1002/cpt.2564

43. Celec P, Vlková B, Lauková L, Bábíčková J, Boor P. Cell-free DNA: the role in pathophysiology and as a biomarker in kidney diseases. *Expert Reviews in Molecular Medicine.* 2018;20:e1. doi:10.1017/erm.2017.12

44. Gauthier VJ, Tyler LN, Mannik M. Blood clearance kinetics and liver uptake of mononucleosomes in mice. *J Immunol.* 1996;156(3):1151-1156. doi:10.4049/jimmunol.156.3.1151

45. Mercier F, Kerioui M, Desmée S, Guedj J, Krieter O, Bruno R. Longitudinal analysis of organ-specific tumor lesion sizes in metastatic colorectal cancer patients receiving first line standard chemotherapy in combination with anti-angiogenic treatment. *J Pharmacokinet Pharmacodyn.* 2020;47(6):613-625. doi:10.1007/s10928-020-09714-z

46. Kerioui M, Desmée S, Mercier F, et al. Assessing the impact of organ-specific lesion dynamics on survival in patients with recurrent urothelial carcinoma treated with atezolizumab or chemotherapy. *ESMO Open.* 2022;7(1):100346. doi:10.1016/j.esmoop.2021.100346

47. Moss J, Magenheim J, Neiman D, et al. Comprehensive human cell-type methylation atlas reveals origins of circulating cell-free DNA in health and disease. *Nat Commun.* 2018;9(1):5068. doi:10.1038/s41467-018-07466-6
48. Mattox AK, Douville C, Wang Y, et al. The origin of highly elevated cell-free DNA in healthy individuals and patients with pancreatic, colorectal, lung, or ovarian cancer. *Cancer Discov*. 2023;13(10):2166-2179. doi:10.1158/2159-8290.CD-21-1252

49. Tsiatis AA, Davidian M. Joint Modeling of Longitudinal and Time-to-Event Data: An Overview. *Statistica Sinica*. 2004;14(3):809-834.

50. Savic RM, Karlsson MO. Importance of Shrinkage in Empirical Bayes Estimates for Diagnostics: Problems and Solutions. *AAPS J*. 2009;11(3):558-569. doi:10.1208/s12248-009-9133-0