Physics-Regularized Multi-Modal Image Assimilation for Brain Tumor Localization

Authors: michal  balcerak
tamaz  amiranashvili
andreas  wagner
jonas  weidner
petr  karnakov
johannes c paetzold
ivan  ezhov
petros  koumoutsakos
benedikt  wiestler
bjoern  menze
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 15:28:57

Abstract: Physical models in the form of partial differential equations serve as important priors for many under-constrained problems. One such application is tumor treatment planning, which relies on accurately estimating the spatial distribution of tumor cells within a patient’s anatomy. While medical imaging can detect the bulk of a tumor, it cannot capture the full extent of its spread, as low-concentration tumor cells often remain undetectable, particularly in glioblastoma, the most common primary brain tumor. Machine learning approaches struggle to estimate the complete tumor cell distribution due to a lack of appropriate training data. Consequently, most existing methods rely on physics-based simulations to generate anatomically and physiologically plausible estimations. However, these approaches face challenges with complex and unknown initial conditions and are constrained by overly rigid physical models. In this work, we introduce a novel method that integrates data-driven and physics-based cost functions, akin to Physics-Informed Neural Networks (PINNs). However, our approach parametrizes the solution directly on a dynamic discrete mesh, allowing for the effective modeling of complex biomechanical behaviors. Specifically, we propose a unique discretization scheme that quantifies how well the learned spatiotemporal distributions of tumor and brain tissues adhere to their respective growth and elasticity equations. This quantification acts as a regularization term, offering greater flexibility and improved integration of patient data compared to existing models. We demonstrate enhanced coverage of tumor recurrence areas using real-world data from a patient cohort, highlighting the potential of our method to improve model-driven treatment planning for glioblastoma in clinical practice.

Main paper: 
 introduction: 
 The management of gliomas, particularly glioblastomas, is highly challenging due to their infiltration beyond the tumor margins detectable in medical imaging. This infiltrative behavior complicates the precise personalization of radiotherapy. Current clinical practice applies uniform radiation to a 1.5 cm margin around the active tumor core visible in Magnetic Resonance Imaging (MRI) scans. However, this approach overlooks critical factors such as heterogeneous infiltration patterns, varying brain tissues, and anatomical barriers, failing to account for the highly individualized dynamics of tumor spread in each patient. Personalizing radiotherapy to account for this spread remains an unmet clinical need in neuro-oncology [1] .The computational modeling of tumor growth has become a promising solution to this challenge [2, 3, 4, 5] , by formalizing the tumor growth process through mathematical models of varying complexity [6] . However, personalizing these models to the limited clinical observational (image) data of a patient remains a formidable challenge. To address it, two main computational approaches have emerged:Hard-Constrained Physics-Based Models: These models use Partial Differential Equations (PDEs) to constrain the solution space, ensuring that output tumor cell distributions adhere strictly to predefined physical laws [7, 8, 9, 10, 11] .However, since PDEs are approximations of the underlying stochastic processes in complex biological environments, they may overconstrain the system. While these models offer structure and insight, they fail to capture the intricate behaviors inherent to biological processes like tumor growth. Additionally, linking image observations to state variables requires careful consideration. As a result, the full integration of data is hindered by these limitations, reducing the effectiveness of physical models in representing empirical data.Data-Driven Models: These models [12, 13] offer great flexibility by leveraging access to imaging datasets. They can directly be learned from -and applied to -the patients' diagnostic scans. However, they often lack rigorous quality control and fail to integrate crucial biological insights, such as tissue-specific infiltration patterns and brain topology. This makes them less reliable for clinical application, especially when there is insufficient data to infer them purely from imaging information.In Search of a Hybrid Approach: Combining physical insights on spatiotemporal distributions of state variables with increased model flexibility could be highly beneficial. One approach to achieve this is by testing the physics residual in selected locations, as demonstrated in Physics-Informed Neural Networks (PINNs) [14, 15, 16, 17, 18] . Another way involves projecting the learned solution into a latent space that adheres to physics constraints [19, 20] . A different framework, Optimizing a Discrete Loss (ODIL) [21] , constructs a residual of PDEs using grid-based discretization. Unlike methods that optimize the weights of Multilayer Perceptrons (MLPs) that implicitly parametrize the unknown fields in the PINNs method, ODIL directly optimizes the discretized unknown fields, which can be faster and more stable due to the local (rather than global, as in MLPs) dependencies of the solution on the learnable weights.Both ODIL [22] and PINNs [23] have recently been employed to model growth equations conditioned on MRI and metabolic imaging of patients. ODIL advances state-of-the-art radiotherapy planning by relaxing the constraints of the growth equation, thereby enhancing the capture of tumor recurrence. A similar trend is observed in computer vision's novel view synthesis, where neural radiance field MLPs [24] are being replaced by discrete representations [25, 26, 27] .However, none of these hybrid approaches have modeled brain-tumor interactions dealing with constraints imposed by local structural anatomy, such as the relative composition of brain tissues and their microstructure, or by the impact of gross anatomical deformations visible in structural scans of large tumors. While biomechanical model extensions can address the latter, incorporating tissue elasticity remains challenging due to the unknown initial anatomy and the computational complexity of the tumor growth dynamics. So far, no hybrid approach free from hard PDE constraints has successfully linked biological processes, such as tumor growth, with models of the biological environment, like tissue elasticity. Only a few studies have used PINNs for biomechanical modeling in synthetic scenarios, demonstrating that the displacement field of materials can be learned from observations [28, 29] .The main contributions of this work are as follows:• We construct discrete physics-based residuals as measures of differences between learned spatiotemporal distributions and those predicted by discretized physical equations governing tumor growth and tissue biomechanics. Using a Lagrangian perspective where particles carry tissue intensities, we project these dynamics onto static Eulerian grids to condition the learnable tumor cell and tissue distributions on available observations. Unlike data-driven models, our hybrid method remains robust due to physics regularization while being more flexible at adapting data than hard-constrained numerical simulations. • Our method allows for the estimation of the initial condition of the biomechanical environment.We demonstrate the ability to learn the complex unknown initial condition, which, in our case, is not just the origin of the pathology but also the initial state of the brain tissue anatomy. • We validate our method * on the downstream task of radiotherapy planning using the largest publicly available dataset known to us, achieving new state-of-the-art performance in capturing tumor recurrence. Figure 1 provides an overview of our method, showing how we condition the learnable full spatial tumor cell distribution on the patient data and regularize it using soft assumptions about the physics of tumor growth, tissue elasticity, and initial conditions for both the anatomy and the tumor. In addition to multiparametric MRI scans (T1Gd, T1, FLAIR), we incorporate Fluoroethyl-L-tyrosine PET (FET-PET), a specialized form of Positron Emission Tomography (PET), as an additional imaging modality that captures the metabolic activity of tumor cells. 
physical assumptions: 
 Our model integrates physical assumptions to ensure that the learned tumor cell distribution adheres to established biological and mechanical principles. These assumptions are incorporated into the loss function as regularization terms alongside initial condition assumptions and data terms.We consider the problem of learning the tumor cell distribution within the domain Ω × (0, 1], where (0, 1] represents the assumed time of tumor growth, and the unit cube Ω = [0, 1] 3 describes our spatial domain.The objective is to determine a tumor cell density function c(x, t), (x, t) ∈ Ω × (0, 1] approximately satisfying the reaction-diffusion-advection equation [30] :EQUATIONwhere Dc = ∇ • (D(m)∇c) represents the diffusion of the cells with the tissue-dependent diffusion constant D(m), S(c) = ρc(1 -c) describes the cell growth with a growth factor ρ, and v is a velocity field moving the cells in the interior of the domain due to the tissue's elastic properties and the tumor-induced stresses in the tissue.The brain tissue vector m, which represents the percentage-wise concentration of different tissue types, is defined as m(x, t) = [m WM (x, t), m GM (x, t), m CSF (x, t)], where m WM , m GM , and m CSF denote the concentrations of white matter, gray matter, and cerebrospinal fluid, respectively.The brain tissue vector m is governed by:EQUATIONwhere ⊗ denotes the outer product. We assume Neumann boundary conditions for both c and m. For c, the no-flux condition applies at the borders of diffusive tissues, specifically the combined white and gray matter. For each tissue component in m, the no-flux boundary is assumed to be the constant brain boundary visible in the MRI scans.The tumor induced stresses lead to a displacement field u, such that the velocity is given by v = ∂ t u.We assume quasi-static mechanical equilibrium due to the slow tumor growth rate relative to tissue mechanical responses, leading to:EQUATIONwhere γ regulates the impact of the tumor on the displacement and is a learnable patient-specific parameter. We apply the Neo-Hookean model [31] to represent hyperelastic behavior in biological tissues, where the stress σ is related to the deformation gradient F ij = δ ij + ∂ui ∂xj with the Lamé parameters λ and µ (material properties in Appendix G), averaged across different tissue types to accommodate the heterogeneity:EQUATIONwith J = det(F). The averaged Lamé parameters, λ and μ, are computed based on the proportional contributions of constituent tissues, ensuring the stress tensor accurately reflects the composite mechanical properties of the mixed tissues. This methodology facilitates detailed simulations adhering to the elasticity laws relevant to heterogeneous biological environments.The set of learnable parameters θ Dynamics describing the tumor dynamics, is thus given by θ Dynamics = {D GM , R, ρ, γ}, where D GM and D WM = RD GM are the diffusivities in pure gray and white matter regions, which are averaged to the effective diffusivity with a weighted average based on the tissue proportions.We use two different approaches to discretize the elasticity equations and the reaction-diffusionadvection equations. For the former, we use a time-dependent grid in space, described by the N x points (p n j ) 1≤j≤N X at time n∆t, 0 ≤ n ≤ N t , which initially form a uniform grid on Ω. For the latter, we further partition the domain Ω at the nth point in time intoN x cells Ω n i , 1 ≤ i ≤ N x , such that Ω = Nx i=1 Ω n i .Note that the nodes defining the cell boundaries are carried by the particles p n j , whose movement, governed by tissue elasticity, gradually deforms the cells Ω n i from their initial rectangular shape.The displacement field u n is computed based on the initial (t = 0) and current (t = n∆t) positions of particles p n within the tissue by u n = p np 0 . Using central finite differences on the initial uniform reference grid, we can consecutively approximate ∇u n j , F n j , ∇c n j , ∇ • σ n j , and define the tissue residual:EQUATION)As the particles, representing grid points, move with the tissue's deformation, they inherently accommodate the advection. Thus, the observed displacement of tumor cells attributed to advection is absorbed into the changes in particle positions. This leads to a simplification of the tumor equation when solved in the moving reference frame of the particles.We use an implicit Euler scheme combined with a cell-centered finite-volume discretization to discretize the tumor equation. We denote the boundaries between two cells 1 ≤ i, j ≤ N x by Γ n ij = Γ n ji = ∂Ω n i ∩ ∂Ω n j . Further, the index set of all neighboring cells of Ω n i is given by N i = {j : ∂Ω j ∩ ∂Ω i ̸ = ∅ ∧ j ̸ = i} and assumed to stay independent of n, due to the extend of our deformations (see limitations in Section 7) . For a function c n i at time n∆t at the ith cell, we discretize the diffusion D and reaction S operators by:D[c n i , D n i ] = j∈Ni |Γ ij |D n ij c n j -c n i ∥x n j -x n i ∥, and S[EQUATIONwhere D n ij is the harmonic mean of D n i at the cell boundary, |Γ n ij | denotes the interface area, and |Ω n i | denotes the cell volume. We thus obtain the tumor loss by: EQUATIONWhile a Lagrangian frame representation from the point of view of particles is beneficial for accurately modeling advection, an Eulerian frame is essential for rendering the advected particle states into the image space. Denoting the Eulerian and Lagrangian frames as G and P respectively, an Eulerian field F G and a Lagrangian field F P are related through:EQUATIONwhere i indexes grid nodes, j indexes particles, and w ij represents the weight of the trilinear shape function defined on node i and evaluated at the location of particle j.Question: Does the paper discuss the limitations of the work performed by the authors? Answer: [Yes] Justification: Yes, we have addressed the multi-resolution grid, which operates under the assumption that the problem at hand can be reliably approximated at lower-than-target resolutions.Guidelines:• The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper. • The authors are encouraged to create a separate "Limitations" section in their paper.• The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be. • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations. 
initial assumptions: 
 We model the initial conditions of the tumor and brain tissues by incorporating assumptions into our loss functions, ensuring the model aligns with the expected biological and anatomical state of the brain before and at the onset of tumor growth.We assume that the initial tumor distribution can be represented by a small Gaussian blob. To enforce this assumption, we penalize the difference between F G (u) 0 and a Gaussian function centered at x 0 ∈ Ω. The center of this Gaussian blob, x 0 , is part of the learnable parameters θ Initial = {x 0 , m 0 }. The initial condition loss for the tumor cells is defined as:EQUATIONwhere D 1 = 0.5, D 2 = 0.02 are constants and x i is the position of the center of Ω 0 i .The deformations are clearly visible in the MRI images and are an integral part of understanding growth dynamics. To obtain a plausible estimate of the unperturbed brain, we need to incorporate an anatomical prior.A healthy brain is roughly symmetric between hemispheres, as justified by anatomical and functional observations. Imaging techniques such as MRI reveal symmetrical tissue structures, and electroencephalography (EEG) recordings show symmetrical patterns of electrical activity in the brain [32, 33] . This assumption applies to key tissue types, including white matter, gray matter, and cerebrospinal fluid. We construct a loss function, denoted as L Healthy Anatomy , to quantify the asymmetry of the learned F G (m 0 ). The implementation details are provided in Appendix E. We do not use the unperturbed anatomy explicitly for inference of the tumor cells, however, we need additional constraints to bound γ from Eq. 3 that couples tumor cells with tissue dynamics. While symmetry of the healthy brain is a coarse assumption, we soften it by quantifying the asymmetry in lower resolutions. Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof? 
imaging data reconstruction: 
 Imaging data refers to information directly observed from the patient's MRI and FET-PET 3D scans [34] , including visible tumor parts and surrounding brain tissues. As shown in Figure 1 , this data includes estimated brain tissue probability maps, visible tumor segmentation, and a metabolic map derived from FET-PET scans, focusing on the enhancing core and edema. MRI scans provide segmentation maps [35] , revealing the visible tumor. Using non-rigid atlas image registration with tumor masking [34] , we estimate tissues both outside and hidden by the tumor, including white matter (WM), gray matter (GM), and cerebrospinal fluid (CSF). These observations offer partial information about both the anatomy and the tumor cell distribution.The imaging function (see Figure 2 ) [22, 36] bridges simulated tumor cell densities with visible tumor segmentation and brain tissues, as well as the metabolic map. It depends on parameters for thresholding the tumor cell distribution to obtain segmentation, denoted as θ Imaging , which must be learned.We construct a loss function that minimizes the difference between model predictions and patient observables. Specifically, L Core matches the visible tumor core, L Edema matches the visible edema, and L Metabolic Map matches the metabolic map profile. Tissue intensities carried by particles are set to reproduce the patient's estimated tissues at t = 1. While these intensities are predetermined, particle positions are weakly constrained, so matching the tissues (WM, GM, CSF) to the observed structures outside the tumor is enforced through L Tissues . The constructed penalty terms align our spatiotemporal distributions of tumor cells and tissues with patient data characteristics. Implementation details of the data losses and the imaging function are provided in Appendix D. 
combined loss function: 
 The loss function L Total is constructed to balance the contributions of physics-based regularization, assumptions about the initial conditions, and the imaging empirical data. This facilitates a comprehensive approach to infer the full spatial tumor cell distribution at the time of the patient's data acquisition (see Figure 3 for an overview of the inference process and Figure 2 for the learning process visualization):EQUATION)All α * parameters were determined based on experiments with synthetic data involving both single focal and local multi-focal tumors. Extensive ablation studies are provided in Appendix C, and the details of the synthetic experiment setup can be found in Appendix F.Question: Does the paper provide open access to the data and code, with sufficient instructions to faithfully reproduce the main experimental results, as described in supplemental material?Answer: [Yes] Justification: Both data and code are available to reproduce the main results.Guidelines:• The answer NA means that paper does not include experiments requiring code. • Providing as much information as possible in supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted. 
validation and ablation study of radiotherapy planning: 
 There is no objective ground truth for pre-operative tumor cell distribution, so methods are evaluated through downstream tasks like defining radiation target volumes. An "accurate" method targets areas that later show progression in follow-up MRI, acknowledging that radiation delays but does not eliminate tumor regrowth. An "inaccurate" method spares these areas or targets regions that remain tumor-free. The current "Standard Plan" defines the target volume using a 1.5 cm isosurface around the MRI-segmented tumor core, while our approach generates isosurfaces from learned tumor cell distributions. Our computational model uses a multi-resolution method [37] with additional, coarser grids to accelerate convergence. We employ four levels of grid refinement, resulting in 152,880,048 unknowns.At the finest grid level of 72 × 72 × 72 × 96, each grid point contains four unknowns: three for particle positions and one for tumor cell density.The learning process includes system identification of unknowns: θ Dynamics , θ Initial , and θ Imaging , which parameterize tumor growth dynamics, tissue elasticity, initial conditions, and imaging characteristics.Using the Adam optimizer, convergence takes around 3 hours on an NVIDIA RTX A6000 GPU.We use a publicly available dataset of 58 patients with preoperative MRI scans and preoperative FET-PET imaging, along with follow-up MRI scans at the time of the first visible tumor recurrence [22] .We would like to emphasize that this is the largest publicly available dataset of this kind, i.e., with multi-modal imaging data and follow-up MRI. Its size is common for medical studies and sufficient to motivate clinical trials (see, for example, [38, 39, 40] on clinical investigations into the use of PET for radiotherapy planning). The problem we address is a critical challenge in cancer treatment research that calls for ML solutions. However, data acquisition is costly and time-consuming, and clinical datasets, such as ours with only a few dozen complete observations, are insufficient to support purely data-driven ML approaches (as visible from results reported in Table 1 ). At the same time, it is this inherent limitation that is encouraging us to incorporate physics-based assumptions and physics-informed regularizations into an ML model to effectively address the driving clinical problem.The follow-up images with visible tumor recurrences are not used by our method during the optimization or hyper-parameter search at any point; therefore, the dataset can be considered 100% unseen. The decision not to use any recurrence data for hyper-parameter tuning is to avoid diluting the already small dataset.Our analysis includes a range of models, each categorized by their distinct capabilities:• Numerical Physics Simulations: Employ methods such as Finite Element Method (FEM) and Finite Difference Method (FDM) to model tumor growth dynamics and account for dynamic tissue behaviors. The dynamic parameters are determined through numerous simulations with varying parameters. The brain tissues' initial conditions are based on average brains.• Data-Driven Neural Networks (Unconstrained): Use Convolutional Neural Networks (CNNs) to directly predict likely recurrence locations, representing a population-based data-driven approach. • Data-Driven Neural Networks (Physics-Constrained): Similar to the above unconstrained approach but limited to finding parameters for physics simulations, incorporating hard physics constraints. • Static Grid Discretization: A previous state-of-the-art method that uses Optimizing a Discrete Loss (ODIL) and static grids to penalize the learned spatiotemporal tumor cell distribution by quantifying the discrepancy with the tumor growth equation, keeping brain anatomy static throughout the entire process. • Standard Plan: Apply uniform safety margins around the tumor core, serving as the baseline for volume determination and reflecting current clinical practice post-tumor resection surgery radiotherapy planning.Specific implementations within these categories are explained in Appendix A.Table 1 , specifically the "Recurrence Coverage [%] (Any)" column, and Figure 5a present results where recurrence is defined by any segmentation visible in the follow-up MRI scan. Our method achieved 74.7% recurrence coverage, surpassing the previous state-of-the-art at 72.9% and the Standard Plan at 70.0%. Using a predefined average brain reduced the performance to 73.4%. The physics-constrained data-driven approach slightly outperforms the Standard Plan at 67.1% and significantly outperforms its unconstrained variant at 59.0%. This highlights the advantage of adding physics constraints. Figure 5a shows the distance between 'Greater' and 'Less' categories increased from 26% to 36%.Table 1 , using the "Recurrence Coverage [%] (Enhancing Core)" column, and Figure 5b show results where recurrence is defined by enhancing core segmentation, per RANO guidelines [41] . Our method leads with 89.9% average recurrence coverage, compared to 89.0% for static grid discretization and 87.3% for the Standard Plan. Numerical physics simulations are more robust at 86.2% compared to 84.3% for the physics-constrained data-driven method. The unconstrained data-driven method remains the lowest at 66.8%. Results are close because the enhancing core often occurs near preoperative MRI segmentations, leading to many methods achieving 100% coverage (seen in 'Equal' in Figure 5b ). However, Figure 5b shows an increase in the preference from 19% to 28%.In both recurrence definitions, mean estimate uncertainty is substantial, as Recurrence Coverage for individual patients can vary significantly. These uncertainties are the standard errors of the mean. Nonetheless, Figure 5 shows that even slight improvements in average Recurrence Coverage consistently enhance radiotherapy outcomes compared to the Standard Plan baseline. As shown in Figure 5 , the gap between the "Greater" and "Less" categories represents the number of patients who benefit from our method compared to the standard plan. We observe a consistent improvement with our method, demonstrating superior performance beyond just higher mean recurrence coverage, which could be skewed by outliers. The results are statistically significant, as shown in Table 3 in Appendix B.Question: Does the paper specify all the training and test details (e.g., data splits, hyperparameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?Answer: [Yes] Justification: Yes, all major aspects are provided.Guidelines:• The answer NA means that the paper does not include experiments.• The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them. • The full details can be provided either with the code, in appendix, or as supplemental material. 
conclusion: 
 In this work, we present the first successful approach to a soft-constrained physics systemidentification problem that combines the complex biological process of tumor growth with the biomechanical environment such as elastic brain tissues. Our method integrates physics-based constraints with multi-modal imaging data to enhance tumor treatment planning for glioblastomas. The method balances adherence to observed patient data and physics-based penalties through a unique discretization scheme, serving as a flexible spatiotemporal regularization term.Additionally, our method provides estimates of the initial condition, i.e., tumor-free, pre-deformation anatomy. Such healthy brain anatomy can be utilized in various studies, e.g., for aligning with post-surgical scans to assess the extent of tumor removal and detect complications [42, 43, 44] .Overall, our results show that this approach outperforms previous state-of-the-art techniques in covering tumor recurrence areas, with improved performance on real-world patient data. Finally, we want to point out that relaxing physical model constraints has broad applicability beyond glioblastoma localization. This approach can extend to other real-world problems where rigid physics-based models are limited but still relevant.Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?Answer: [Yes] Justification: Yes -this information is provided.Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).• The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors). • It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text. 
experiments compute resources: 
 Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?Answer: [Yes]Justification: Yes -this information is provided.Guidelines:• The answer NA means that the paper does not include experiments.• The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage. • The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute. • The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper). 
code of ethics: 
 Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https://neurips.cc/public/EthicsGuidelines? 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed? 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?Answer: [NA]Justification: The method is evaluated on a publicly available dataset and does not inherently require safeguards.Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort.12. Licenses for existing assets Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?Answer: [Yes] Justification: Yes -the dataset authors are credited with a licence.Guidelines:• The answer NA means that the paper does not use existing assets.• The authors should cite the original paper that produced the code package or dataset.• The authors should state which version of the asset is used and, if possible, include a URL. • The name of the license (e.g., CC-BY 4.0) should be included for each asset.• For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided.• If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswithcode.com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a dataset. • For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided. • If this information is not available online, the authors are encouraged to reach out to the asset's creators. 
new assets: 
 Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?Answer: [NA] Justification: No new assets. (other than later provided code)Guidelines:• The answer NA means that the paper does not release new assets.• Researchers should communicate the details of the dataset/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc. • The paper should discuss whether and how consent was obtained from people whose asset is used. • At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file. 
crowdsourcing and research with human subjects: 
 Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?Answer: [NA] Justification: The authors were not involved in data gathering. The patients protocols are explained at the cited source.Guidelines:• The answer NA means that the paper does not involve crowdsourcing nor research with human subjects. • Including this information in the supplemental material is fine, but if the main contribution of the paper involves human subjects, then as much detail as possible should be included in the main paper. • According to the NeurIPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector. 
institutional review board (irb) approvals or equivalent for research with human subjects: 
 Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?Answer: [Yes]Justification: The authors were not involved in data gathering. The patients protocols are explained at the cited source.Guidelines:• The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.• Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper. • We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution. • For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review. 
 
 Appendix: The convergence speed of the method depends on a multiresolution grid, assuming that problems can be represented at lower resolutions. However, this approach may not be directly applicable to PDEs with narrow-phase interfaces. Additionally, while large deformations are allowed, we assume they are not large enough to cause material rupture, which may not be optimally handled in the Lagrangian perspective without modifications to the presented scheme.Here we provide details regarding the specific baseline implementations:Numerical Physics Simulations: We employ the numerical forward simulation scheme from [45] which was implemented on GPU [13] . The simulations follows the equation 1. We initialize the brain tissues with a rigidly registered average brain. For θ Dynamics , x 0 ∈ θ Initial , and θ Imaging , we use the Covariance Matrix Adaptation Evolution Strategy (CMA-ES) method [10, 46] .Data-Driven Neural Networks (Unconstrained): We use nnU-Net [47] which takes brain tissue distribution, FET-PET map, and pre-operative tumor segmentations to directly predict core tumor recurrences. After predicting the recurrences, we construct a distance distribution from the predicted recurrences. This distribution is then thresholded, similar to other methods, to create a binary radiotherapy treatment map of the same volume as the Standard Plan. The network is trained using an 80%/20% train-test split.Data-Driven Neural Networks (Physics-Constrained): We employ a recently published method [13] that uses CNNs to map pre-operative tumor segmentation to θ Dynamics and runs a numerical simulation following equation 1. This method uses an average brain as the initial anatomy. After inference, the tumor cell distribution is non-rigidly registered to the patient's brain anatomy.Static Grid Discretization: This method, as implemented in [22] , uses estimated brain tissues (4) as the initial condition for the anatomy and later keeps the tissues static, effectively setting γ ∈ θ Dynamics in Equation 3 to 0.Standard Plan: This approach uses a 1.5 cm uniform safety margin around the resection cavity and/or remaining tumor on MRI imaging, aligning with both North American and European guidelines [48, 49] .To show how the method performs on tasks where there is ground truth, our method and GliODIL [16] (previous SOTA, Static Grid) are additionally evaluated on synthetic cases generated by numerical physics solvers under various conditions. See Table 2 for the results and Table 3 Additional information regarding the loss coefficients:α 2 is essential for grid stability. Without it, the grid would be unstable, causing folding and rendering the results meaningless. If set to zero, α 2 can be replaced by a Laplacian-based regularization, as visualized in Figure 6b and Figure 6c , though this clearly results in less biologically plausible deformations.α 4 regularizes the initial anatomy to bound γ from Eq 3. It should be kept low due to the natural asymmetry in brains. More details are provided in Appendix E.α 5 , α 6 , α 7 , and α 8 are crucial for grounding PDE solutions in reality (most important, all other losses can be seen simply as regularisation of these components). Imaging data reconstruction loss cannot be calibrated using synthetic studies as all the other alphas. Instead of using recurrence data to calibrate them, we set them such that they provide a balanced contribution to the total loss α 3 should be kept low since the initial tumor shape is unknown. However, decreasing it too much causes a local minimum with no tumor.α 1 serves as the main regularization term. It should not be → ∞ to avoid solutions not flexible enough to properly accommodate patient data (see Numerical Physics Simulations results in Table 1 ; this setup performs poorly), nor should it be → 0 to prevent overfitting to imaging data alone, effectively collapsing to standard plan level performance.Removing elements of the imaging data reconstruction loss introduces ill-posedness and degradation of synthetic results (see Table 4 ). Removing tumor growth physics regularizations effectively reduces our solution to performance close to the standard plan. Removing particle dynamics (static grid) reduces our model to [16] Static Grid Discretization ("GliODIL"). The core of the method's alignment with the data is encapsulated by the loss function components, which associate the simulated outputs with key imaging traits such as the tumor core, surrounding edema, metabolic activity detected through FET-PET imaging, and visible brain tissues. This loss function depends on two learnable parameters θ Imaging = {θ down , θ up } and is expressed as:EQUATIONThe loss function is composed of individual terms corresponding to distinct anatomical features:• L Core relates tumor cell concentrations above the threshold θ up to the tumor core region.• L Edema delineates the edema area surrounding the tumor, regulated by the lower and upper thresholds θ down and θ up .• L Metabolic Map assesses the metabolic activity as indicated by FET-PET signals within the edema and core regions, using a simple correlation metric, L Tissues connects visible brain tissues with those inferred by the method.These adaptive parameters {θ down , θ up } enable the model to accommodate variations in MRI/FET-PET imaging contrasts and noise levels.We adopt sigmoid functions to portray the gradational transitions observed at tumor region margins. The sigmoid, σ(x), is specified as:EQUATIONHere, β modulates the steepness of the transition and is set to β = 50. For the tumor core:EQUATIONfor the edema:EQUATIONwhere α offsets the thresholds and is set to α = 0.05.The metabolic activity within the tumor is evaluated by the loss term L PET , which measures the correlation between the simulated metabolic signal and actual FET-PET scan observations:EQUATIONHere, corr(c Nt , p PET ) represents the Pearson correlation coefficient between the predicted tumor cell densities c Nt and the restricted FET-PET signal p PET , across all voxels within within the edema and enhancing core regions, providing a direct measure of how well the model predictions align with the observed metabolic profiles.Our model employs a multi-resolution analysis to assess brain symmetry, particularly useful since the brain's symmetry is more evident at coarser resolutions. This approach is applied to key tissue types: white matter (WM), gray matter (GM), and cerebrospinal fluid (CSF), each analyzed at full, 2x, and 4x downsampled levels to accommodate both the brain's structural nuances and computational efficiency.The process for each tissue type involves converting tissue data into particles, then calculating symmetry loss at each resolution. This calculation entails:1. Downsample the tissue representation, using average pooling 2. Divide the tissue's downsampled representation along its height, reflecting the brain's hemispherical division. 3. Mirror one hemisphere and compare it to the other, quantifying symmetry by calculating the mean absolute difference between them.The total symmetry loss accumulates across all tissues and resolutions to measure deviation from an ideal symmetric state:EQUATIONThe multi-resolution approach not only addresses the limitations of assuming brain symmetry at detailed levels but also streamlines computational efforts. The smoother landscape at downsampled resolutions is especially favorable for particle movement through gradient descent.To calibrate the loss function weights, we generated a synthetic dataset for single focal and localized multi-focal tumors by solving PDEs using a finite difference method numerical solver [13] . An average brain model was used to represent the spatial distribution of brain tissues. Tumor growth model parameters, elasticity, imaging model parameters, and 1 or 3 focal locations were varied using uniform random distributions.The metric used was the RMSE between the learned and ground truth tumor cell distributions. In total, 100 synthetic patients were generated. Table 5 provides the details of the parameter ranges used for generating the synthetic single focal and multi-focal tumor datasets. , y 2 0 , z 2 0 ) (x 1 0 , y 1 0 , z 1 0 ) ± 9.6 Tumor 3 center (x 3 0 , y 3 0 , z 3 0 ) (x 2 0 , y 2 0 , z 2 0 ) ± 9.6The table below summarizes the Young's modulus and Poisson's ratio for different tissue types used in the hyperelastic modeling of biological tissues. Guidelines:• The answer NA means that the abstract and introduction do not include the claims made in the paper. • The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers. • The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings. • It is fine to include aspirational goals as motivation as long as it is clear that these goals are not attained by the paper.Justification: We are here to help, no cause harm.Guidelines:• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction).Justification: There are no clear societal impacts, except for the potential for improved healthcare, which is difficult to assess without clinical trials. The work also has the potential to translate to other areas.Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.• Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 


Decision: This work presents a new approach to integrate physics-based tumor growth constraint with multi-modal imaging data to predict tumor cell distribution and thus enhance tumor treatment planning for glioblastomas. The core of the methods includes a discrete physics residual and initial assumptions encoding initial tumor distribution and symmetrical pattern between hemispheres. Experiments were performed on a small real dataset, and the success was measured by the recurrence coverage defined as the percentage of tumor (based on follow-up MRI) within the model-predicted volume. Experimental results demonstrate improved recurrence coverage of the presented method compared to fully data-driven or fully physics-based methods.

The presented work is of important potential in bringing physics-informed learning into biomedical tasks and demonstrates the ability to integrate physics-based constraints and multi-modal imaging data. The proposed method of discretizing the physics residuals is novel and potentially applicable to related problems and of interest to the research community. The paper is clearly written and the technical details are sound and sufficient.

The reviewers have some concerns: the experimental evaluation of the presented work was performed on a very small dataset with unclear training-test split. This raises some question about the general conclusion that can be drawn from the results obtained on limited samples. Given the reported standard deviation, the obtained margin of improvement appears to be marginal without statistical significance. The objective of the presented model includes a large number of terms. The sensitivity of the model performance to the include/exclusion of these different terms and their hyperparameters deserves substantial analyses that are missing in the current paper. Further concerns are about the evaluation: the main weakness is the lack of ground truth data to train and evaluate the model performance. The evaluation is performed not on the original task of brain tumor detection but rather on a downstream task of “recurrence coverage”, which measures the percentage of the tumor detected in follow-up MRIs (rather than the original). Related to above, the authors only use one metric to evaluate their model, which may be insufficient, especially since there is no ground truth labels. 

The rebuttal addressed these concerns, and the paper has been improve accordingly. Most reviewers are convinced by the rebuttal.  The work addresses an important research problem and shows great promise, it has the potential to be highly impactful.

Reviews: review: paper_summary: paper_summary: This paper presents a method for brain tumor modeling through a joint data-driven and physics-based approach. The proposed approach leverages a multi-task loss that incorporates physical assumptions and prior healthy anatomy to more faithfully model the tumor cell distribution. This method outperforms existing computational approaches and the clinical standard with respect to recurrence coverage
strength_weakness: strength_weakness: - I am slightly confused by the recurrence coverage metric. Has this metric been used in prior work? Why not a well-established overlap metric such as intersection over union, for example?
questions: questions: - According to line 207, recurrence coverage is the “percentage of the tumor segmentation… that is encompassed within the radiotherapy target volume defined by each model.” If this is true, then wouldn’t a predicted radiotherapy target volume consisting of the entire scan obtain 100% recurrence coverage? Unless I am misunderstanding, an existing overlap metric (which penalizes such “overprediction”) might be more appropriate.
- Line 5: “Deep-learning based” -> “Deep learning-based”
- Line 234: Missing parenthetical “(Physics-Constrained:”
limitations: limitations: Limitations are briefly addressed in Section 7.
paper_summary: paper_summary: This paper presents a new approach to integrate physics-based tumor growth constraint with multi-modal imaging data to predict tumor cell distribution and thus enhance tumor treatment planning for glioblastomas. The core of the methods includes a discrete physics residual and initial assumptions encoding initial tumor distribution and symmetrical pattern between hemispheres. Experiments were performed on a small real dataset, and the success was measured by the “recurrence coverage” defined as the percentage of tumor (based on follow-up MRI) within the model-predicted volume. Experimental results demonstrate improved recurrence coverage of the presented method compared to fully data-driven or fully physics-based methods.
strength_weakness: strength_weakness: The experimental evaluation of the presented work was performed on a very small dataset with unclear training-test split. This raises some question about the general conclusion that can be drawn from the results obtained on limited samples.

Given the reported standard deviation (eg. Table 1), the obtained margin of improvement appears to be marginal without statistical significance. 

The writing of the paper lacks clarify at places. For instance, the results presented in Figure 5 was unclear. What are the “greater”, “less”, and “equal” categories, and why a bigger gap between the “greater” and “less” categories indicate a favorable performance?

The objective of the presented model includes a large number of terms. The sensitivity of the model performance to the include/exclusion of these different terms and their hyperparameters deserves substantial analyses that are missing in the current paper.
questions: questions: Clarifications on the data split used to derive the reported results will be appreciated. 

The statistical significance of the reported performance (Table 1) will be appreciated. 

Clarifications on the results related to Fig 5 (see questions above) will be appreciated. 

How sensitive are the presented methods to the different terms in the loss functions, and which are the most important terms?

How is asymmetry measured as an optimization objective?
limitations: limitations: The authors presented some discussion about the limitations of the current work. Adding discussion about limitations related to the limited sample size for the experiments would be appreciated.
paper_summary: paper_summary: The authors propose a method using soft physics-based regularization to predict the distribution of brain tumor. The main contribution is a novel discretization scheme of physics equations to model brain tumor.
strength_weakness: strength_weakness: 1. The main weakness is the lack of ground truth data to train and evaluate the model performance. The evaluation is performed not on the original task of brain tumor detection but rather on a downstream task of “recurrence coverage”, which measures the percentage of the tumor detected in follow-up MRIs (rather than the original). 

2. Related to above, the authors only use one metric to evaluate their model, which may be insufficient, especially since there is no ground truth labels
questions: questions: 1. Is there any bias in the evaluation metrics since there are no ground truth labels and authors use recurrence coverage?
limitations: limitations: 1. The authors should discuss the limitations of their methods in more details, especially related to the lack of ground truth labels to train the model
paper_summary: paper_summary: The authors proposed a discrete physic-based loss function which can be minimized to infer the spatial distribution of tumor cell distribution, which can not be assessed by imaging data only. The loss maps imaging data to a tumor tissue biomechanical model to solve inverse problems and to better understand individual tumor progression processes. One specific promise of the model is that it can learn the unknown initial condition of the brain tissue before the pathology. Their model is validated by achieving results comparable to state-of-the-art results on the task of tumor recurrence prediction from preoperative brain MRI.
strength_weakness: strength_weakness: In the reference [16], presented as the previous state-of-the-art method for recurrence coverage prediction, it seems that they obtained better results than the results displayed in table 1. More precisely, they obtained a recurrence coverage on all tumor areas of 76.09% which is better than the method presented in the submitted paper. It appears to me that the models were validated on the same dataset, but perhaps I missed something. Could the authors clarify this point?

Following the previous point, it seems the paper lacks discussion about the use of the “Healthy anatomy” loss term. Compared to the results displayed in [16], it seems that there are no significant improvements using this new loss. How can we assert the correctness of the initial healthy brain tissue reconstruction? The authors claim that this initial state could be used to better identify post-surgical complications, but no clear method is described.
questions: questions: - Could the authors highlight the differences between [16] and their paper? Does addition of the tissue deformation in the physical modeling really improve the results?
- Could the authors better explain the interest of the initial condition reconstruction?
- Is the data-driven neural network trained on only 20% of the 58 data from radiotherapy planning dataset? Is it enough to train such a model? Could we obtain better results by increasing the training dataset?
- How are inferred the parameters theta_up and theta_down? 
- In the appendix, in section B, it seems that c_i should be written instead of u_i (equation 14).
- While it seems that visually the brain tissues exhibit symmetry, it should not be completely symmetric as we can observe on Figure 3 on the average brain template. It seems that “Healthy anatomy” term could lead to overly symmetric structures (see Figure 3 on the learned initial condition between the right and left hemisphere). It appears to me that this could be a limitation of the reconstruction method. Could the authors elaborate on this? The value of the paper could be increased with some analysis/discussion of this reconstruction.
limitations: limitations: The authors incorporated a “limitations” paragraph to address the convergence seed of their method and their tissue deformation model. However, the paper lacks some discussion in regard to their contribution compared to reference [16]. To which extent does the tissue deformation model improve the results? How can we assert the correctness of the initial condition reconstruction?

Comments: title: HP search, a general scheme
comment: We appreciate the reviewer's valuable feedback and are grateful for updating the score. Please feel free to reach out if any further questions arise.
Thank you for your valuable feedback and for updating the score. The balance between the importance of PDEs and data is indeed a broader challenge in the physics-informed ML field. The reviewer rightly highlights a crucial issue in this area. For such an ML problem, it is essential to use a diverse set of complementary assumptions and regularizations because, without them, the problem becomes ill-posed. All the terms we included are well-grounded with clear physical motivation. 

Regarding sensitivity, the solution should not be overly sensitive to individual changes. We refer to Fig. 2a (and to some extent, Fig. 2e) in [16] (Static Grid Discretization method), where the authors observed a 20% difference in DICE scores by varying $\lambda_{PDE}$ across five orders of magnitude (which corresponds to $\alpha_1$, in our case), which can be seen as low sensitivity to small changes. 

When it comes to actual values for hyperparameter selection, in our work, we opted for a general approach using calibration on synthetic data. We began by determining the order of magnitude for the relevant hyperparameters (e.g., ensuring no tissue folding and preventing the initial tumor cells from vanishing) through a grid search, followed by further fine-tuning using the same grid search method and the RMSE as the metric (as in Table 2 of the one-page PDF). While we do not claim that the parameters we identified are necessarily optimal, they are effective enough to outperform all other baselines without raising concerns about an overly complex hyperparameter optimization scheme that might not generalize to other datasets. For datasets involving pre-operative tumor images, our method should work out-of-the-box, as it relies solely on adherence to the pre-determined growth equations and tissue elasticity, without needing specific information from our dataset. We will extend Appendix D to provide more detailed information on the hyperparameter selection procedure involving grid search and the individual physical motivations behind each parameter to assist researchers in applying our framework to new domains/equations. We are always happy to consult further on this.

As a final point, we would like to emphasize that almost all physics-informed ML frameworks rely on implicitly optimizing the unknown fields using CNNs or MLPs, which introduces numerous additional hyperparameters for these architectures. In contrast, our approach only requires setting the resolution aligned with the data's resolution and the stability conditions of the corresponding growth PDEs, thereby eliminating the need to include these parameters in the hyperparameter search.
I'd like to thank the authors for the effort put in the rebuttals, which has helped addressed some of my previous concerns especially in adding significance test to the experimental results. I do acknowledge the value of the study especially the use of physiological terms to enable optimization given small data. I do however believe that the sensitivity of the method performance to the various hyperparameters used in the loss function deserves a much more comprehensive study -- the rebuttal did a good job demonstrating how the absence/presence of these individual terms would affect the model performance, as well as their rough range of values (large or small). What need to be demonstrate however are the best practice for tuning these individual terms, and how their values affect the margins of improvement the method is currently demonstrating, i.e, how sensitive is the demonstrate margin of improvement to the tuning of the many hyperparameters in the loss, and what would be the strategy to obtain such optimal hyperparameters in real world settings when given a new dataset. I did raise my rating although I think more work is needed to clarify this issue considering the many hyperparameters and their effect on the margin of improvements demonstrated.
Thanks for the detailed answer.

I updated my rating accordingly. I believe the paper is stronger with those changes.
We appreciate the reviewer's valuable feedback and are grateful for raising the score from 5 to 6. Please feel free to reach out if any further questions arise.
We want to thank the reviewer again for the valuable suggestions. 
If we have now addressed all your concerns and questions, we would appreciate if you could consider updating your scores to reflect your post-rebuttal opinion.
Please kindly let us know if you have any follow-up questions or areas needing further clarification.
Dear Reviewer S9Sg,

We thank you again for taking the time to read our manuscript and providing valuable feedback. We hope the rebuttal and additional experiments we provided were helpful. We hope the situation regarding [16] and the initial condition importance is now clear.
Please let us know if there are still pressing concerns.

Thank you, 
The Authors
Dear Reviewer YNED,

We thank you again for taking the time to read our manuscript and providing valuable feedback.
We hope our rebuttal and the additional experiments were helpful in clarifying the data splits, demonstrating the statistical significance of our findings, and highlighting the importance of each loss term.
Please let us know if there are still pressing concerns.

Thank you, 
The Authors
I thank the authors for carefully addressing the reviewers' comments and improving the paper accordingly. I do not have further questions or suggestions on the paper. 

I updated my recommendation from "borderline accept" to "weak accept" as I believe the paper is stronger after the revisions.
I acknowledge that I have read the authors' rebuttal and appreciate their hard work. In particular, thank you for clarifying the recurrence coverage metric and including other overlap metrics.