A Brief Introduction to TOPAS

Improvements to radiotherapy and imaging can be achieved by understanding how subatomic particles travel through apparatus and tissue. The most precise calculations of such phenomena follow the Monte Carlo (MC) method. The TOol for PArticle Simulation (TOPAS) software project, launched in 2009, has helped to greatly improve the use of Monte Carlo simulations for cancer research and treatment (1). Requiring no programming knowledge by its users, TOPAS provides a flexible framework to design simulations for radiation therapy. It enables both clinical applications (e.g., high-precision patient dose calculation) and research (e.g., four-dimensional time-of-flight simulations for detector developments), while its design promotes inter-institutional collaboration (2–5).

One of the main reasons for the success of TOPAS is the parameter control system at its core. Parameters specified in one or more text files define the properties of the simulation. For instance, when used for proton dose calculations, they define the extent of the simulated geometry along each cartesian coordinate, the number of voxels in each direction of the patient computed tomography (CT) scan, the beam angle and the scanning pattern of the beam. The hierarchical nature of these text files allows for sharing and adaptation across radiotherapy centers that use similar equipment (treatment machines, detectors, etc.). The parameter control system syntax is described in detail on the “Welcome to the TOPAS documentation!” ( https://topas.readthedocs.io/en/ latest/) website and by Perl et al. (1).

TOPAS has been expanded from its initial focus on medical physics to also cover radiation biology, offering scoring options that include linear energy transfer (LET), multiple relative biological effectiveness (RBE) models for proton therapy and organ-based outcome models [tumor control probability (TCP) and normal tissue complication probability (NTCP)] (6, 7). TOPAS is now part of the U.S. National Cancer Institute's (NCI) Informatics Technology for Cancer Research (ITCR) Initiative.

An extension system was developed for TOPAS that allows customization for users who need additional features specific to their own application, such as custom scorers or geometries. The extension manager allows users to add new features by writing a short section of simple, very specific C++ classes (a building block of the code) based on templates and C++ helper functions provided within TOPAS. This C++ code can then be integrated into the main TOPAS executable by running a single CMake ( https://cmake.org/) command.

Documentation of the latest version of the TOPAS system and the available parameters can be found at “Welcome to the TOPAS documentation!” ( https://topas.readthedocs.io/en/latcst/).

Motivation for TOPAS-nBio

TOPAS has been successfully applied to research in radiation therapy physics and macroscopic organ or cellular biology. However, more fundamental research is needed to understand the underlying mechanisms of radiation action, describe effects of oxygenation, intracellular signaling, drug-induced radiation sensitization or resistance and many other effects (8–10). Variations at the (sub-)cellular (cellular and sub-cellular) level, both for tumors and surrounding normal tissues, need to be considered. Such research is ideally supported by detailed in silico simulations at the sub-cellular level (11).

The goal of this work was to lay the foundation to a deeper understanding of the biological effects of radiation in order to facilitate new research at the boundary between physics, chemistry and biology. By providing detailed physics and chemistry simulations in combination with detailed representations of biological systems, such as cells and their nuclei, we aim to promote a mechanistic description connecting sub-cellular energy deposition phenomena to observable biological outcomes. We have thus developed TOPAS-nBio, an extension to TOPAS specifically aimed at the simulation of radiobiological experiments by modelling detailed biological effects at the nanometer scale. By taking advantage of the simplicity and reliability engineered into TOPAS while providing nanometer-scale Monte Carlo simulations, we have made complex code accessible to researchers who may consider using Monte Carlo simulations to improve the physical, chemical and biological description of their experimental design or data analysis.

TOPAS-nBio Extensions

The TOPAS-nBio extensions provide options for sub-cellular geometries, scoring, physics and chemistry. While most users will just need to adjust the values of the parameters used by the new classes, the classes can be modified with minimal coding requirements to adjust each aspect of a user's simulation. Due to the modular nature of the extensions, users can select to only install the features they need for their simulations by downloading the necessary files and adding them to their TOPAS executable.

TOPAS-nBio facilitates and extends the use and configuration of the physical and chemical processes provided by Geant4-DNA (15–18). The physical processes of Geant4-DNA, originally intended for radiation transport in liquid water, have recently included cross sections for DNA constituents (19), due to the availability of elastic and inelastic cross sections for these materials, which are also available in TOPAS-nBio.

The Geant4 toolkit provides users with basic 3D geometric shapes (solids), which include volumes such as boxes, ellipsoids, cylinders and spheres. All geometries in TOPAS-nBio are made up of either a single Geant4 solid, computer-aided-designed (CAD) solids or a combination of two or more of these solids. In some cases, unions of these solids are used (e.g., to create the double-helix DNA backbone from the union of spheres). The Geant4-DNA user community and others working on track-structure codes also work on new DNA or cell geometries and new DNA repair models. We actively collaborate with some of these groups to include the latest developments in our simulation framework. Several of the features presented here, e.g., interfaces to the DNAFabric code (20–22) and DNA repair models (23–26), were developed as part of such collaborations.

Physics and Regions

For Monte Carlo track-structure simulations, the transport of charged particles and their interactions are performed in a step-by-step fashion by using the physical processes, their valid energy ranges and associated models provided by Geant4-DNA, described elsewhere (15, 17, 19, 27–30). An overview of the performance using different sets of the most recent physics models available in Geant4-DNA is given by Incerti et al. (15). Physical processes are available through Geant4 constructors. A constructor is a C++ class which implements a physics list, composed of a specific set of physical processes and associated models, configured according to the energy limits of the models and the particle types associated with the processes. Geant4-DNA currently offers seven physics list constructor options: the default constructor, plus options 1 through 6 that offer different scattering models. Geant4-DNA provides cross sections for electrons, photons, protons and a few selected ions (H through O, Si and Fe). Photon cross sections are based on the Livermore models in Geant4. Currently, Geant4-DNA physics cross sections are available for water and DNA-related materials, i.e., tetrahydrofuran (THF), trimethylphosphate (TMP), and pyrimidine (PY) and purine (PU), which serve as precursors for the DNA (and RNA) (19). The only other material that has dedicated low-energy cross sections in Geant4-DNA is silicon (31, 32) using the “MicroElec” package (33, 34). Gold cross sections are expected to be available in Geant4-DNA soon (35, 36).

TOPAS-nBio gives access to the physics list defined in constructors through parameter settings, while providing flexibility to control the model type involved in each process. Using this approach, the energy cut for applying electron capture or electron solvation is automatically readjusted according to the lower energy limit of the physical models. This makes it possible to combine the elastic models from the CPA100 implementation (28) available in constructor G4EmDNAPhysics_opt6 with the inelastic models from the Emfietzoglou-based implementation (30) available in G4EmDNAPhysics_opt4.

Regions were originally developed in Geant4 to allow the use of different physics settings (e.g., production cuts of secondary particles) in different parts of the simulated World, i.e., a volume encompassing the entire simulated geometry. This has been exploited in Geant4-DNA to use different physical models in different geometry components (37), such as limiting tracking of detailed particle interactions to specific components, in order to speed up the simulation. This feature has been implemented in TOPAS-nBio. An example of parameter settings used to simulate two regions with different physics settings, together with the resulting simulation of a gold nanoparticle (GNP) in water, is shown in Fig. 2.

FIG. 2.

The use of different models in different regions. Panel A: Parameter file setting physics lists for different regions to combine macroscopic and nanometer scale simulations with several parameters that should be considered. Panel B: Protons (1 Me V) traversing a gold nanoparticle (GNP) surrounded by liquid water. In water, the transport is handled by Monte Carlo track-structure simulations, whereas in the GNP the transport is handled by condensed-history Monte Carlo (Geant4 Livermore EM physics processes). Proton tracks are shown with blue lines and electron tracks are shown with red lines connected with yellow points that represent inelastic and elastic interactions.

Simulation Settings and Scenarios

The main difference when setting up simulations using TOPAS-nBio, compared to TOPAS, is the setting of specialized physics lists and chemistry lists to define the basic interaction properties for the simulations. In addition, some settings are optimized by default for all TOPAS-nBio simulations to improve the overall performance. For example, the Geant4 geometric tolerance is adjusted depending on the size of the World or can be adjusted by the user. The tolerance defines an envelope around a surface (similar to an uncertainty band) at which Geant4 considers a point to be on the surface. This is a critical parameter at nanometer scales.

The following are typical, distinct simulation scenarios:

Time Dependence

The TOPAS time feature system (4) is used to handle time dependencies within the TOPAS-nBio extension. All time features available in the standard TOPAS are available; that is, beam parameters can vary over time, geometries can move over time, scorers can be switched on or off for certain times, and so on. For example, sources can be defined to have different dose rates, influencing the time between potential cellular damage, which can be important information influencing repair kinetics in the biological effect models linked to TOPAS-nBio, e.g., by including the damage time in the SDD format (see sub-section, “The standard for DNA damage scorer”).

Code Maintenance

TOPAS-nBio was developed by a core team. In addition, alpha users (e.g., University of Manchester) contribute code to the project. Code reliability is fundamental in the development of TOPAS. TOPAS-nBio uses the same system of rigorous development procedures as TOPAS. Each time code is committed to the GitHub repository, a series of instant tests are initiated, including one that checks for a successful build, several “sanity checks” that guarantee the committed code does not fundamentally break the code and a simple regression test, ensuring consistency at several key aspects. In addition, before a new version release, each committed code is checked to conform with internal coding guidelines and a more extensive regression test is initiated. This also occurs every time TOPAS is linked to a new Geant4 release. TOPAS-nBio itself is an extension to TOPAS, thus most TOPAS modules are independent of TOPAS-nBio. In addition, similar to TOPAS releases, TOPAS-nBio releases will check each module for coding consistency and initiate a set of sanity checks and regression tests before code is added to the repository.

Geometries

TOPAS-nBio offers the user a catalogue of geometries ranging from the micrometer scale (e.g., cells) down to the nanometer scale (e.g., DNA strands). Each geometry has its own set of parameters that can be controlled by the user in the parameter file. For example, the cell geometry has parameters for the size of the cell and also the option to include organelles such as a nucleus. Advanced users also have the option of editing or writing their own geometry extension classes to define new unique geometries. A detailed description of the TOPAS-nBio geometries can be found in (40).

Possible targets in radiobiological studies may include the whole cell and nucleus or a specific molecule. A collection of geometries is included to cover a range of anticipated simulations, including several cell morphologies, organelles, nuclei models, various DNA models and proteins. An example parameter file setting up a chromatin fiber using the TsFiber geometry, described in detail by Henthorn et al. (26), is shown in Fig. 3 and a representation of the geometry is shown in Fig. 4C.

FIG. 3.

Parameter file to set up a chromatin fiber defined by a class called TsFiber.

FIG. 4.

Geometries available in TOPAS-nBio. Panel A: An ellipsoid cell shown with a nucleus (blue) and mitochondria (red), a spherical and a fibroblast cell with nucleus and mitochondria. Also shown is a hippocampal neuron with a soma (red) and dendrites (black and blue). Panel B: Three full nucleus models, one based on the Geant4-DNA example (left side) and two different fractal models (center and right side). Panel C: A chromatin fiber consisting of nucleosomes each composed of histone proteins (blue) wrapped by two turns of a double helix DNA (green and red). Panel D: A circular plasmid consisting of 100 basepairs. Panel E: RNA strand recreated using the TOPAS-nBio interface to the protein database. Panel F: A lipid (membrane) layer.

Available cell geometries include simple cells such as spherical or ellipsoid cells, as well as highly specialized cells such as neurons (Fig. 4A). Neuron or glial cell geometries are modeled via the TOPAS-nBio interface to the NeuroMorpho database,  www.neuromorpho.org (41, 42). Since radiation-induced damage to DNA is an important measurement in radiobiology studies, several published DNA models are available, including three full nuclear DNA models incorporating the full hierarchal folding scheme of the DNA in the nucleus (Fig. 4B).

Variance Reduction with Flagged Uniform Particle Splitting

A method to reduce simulation time, in addition to optimization of regions as discussed in the sub-section, “Physics and regions,” is variance reduction. Population control methods such as particle splitting and cross-section enhancement are examples of variance reduction techniques that have been successfully applied to condensed-history simulations with an impressive reduction in computation time [e.g., (43, 44) and papers therein]. In such cases, the quantity of interest is averaged over multiple histories where the track structure of each history is not considered. For track-structure simulations, however, the stochastic nature of each single history is paramount and the artificial generation of new tracks within the same history potentially biases the results. For example, the use of particle splitting will result in an overestimation of the cluster sizes of ionizations in small volumes and the frequency of DSBs. To exploit the large benefits of variance reduction without bias, we implemented a flagged uniform particle split in TOPAS-nBio, specifically designed for track-structure simulations. This technique performs uniform splitting (45) to secondary electrons produced in ionization events at strategically located regions within the geometry and assigns a unique flag number, which is inherited by their progeny (46). The flag permits reclassification of each split event as if they were produced by independent histories. This method reduces the variance by improving the statistics of secondary electrons, while keeping the time increase small compared to the generation of additional particles, by only producing them in strategically selected regions. As a result, improvements of the computational efficiency up to a factor of approximately 65 can be achieved without loss of accuracy (46).

Chemistry

Step-by-step method. Geant4-DNA is capable of simulating the physical, pre-chemical and chemical stages of the interaction of radiation with biological tissue in a single application (47, 48). Configuration of all the simulation parameters can be cumbersome due to the large number of variables involved: dissociation channels for water and their probabilities, reaction rates and diffusion coefficients for each and every reaction, along with configuration of the simulation setup for validation against measured data. TOPAS-nBio provides functionality that facilitates the configuration of simulations that include chemistry (38). TOPAS-nBio further allows expansion of the Geant4-DNA default reaction database and type of chemical species. In addition, TOPAS-nBio provides a set of specialized scorers for the spatiotemporal information of chemical species and their yield through a G-value scorer. The chemical parameters were updated from those used by Geant4-DNA and verified through a comparison with experimental data, obtaining satisfactory agreement, as shown in Fig. 5 and (38).

FIG. 5.

TOPAS-nBio calculated G values as a function of LET_100eV for mono-energetic electrons (e^-), protons (p) and alpha (α) particles. Point-to-point differences with experimental data are shown at the bottom of each panel with error bars (1 standard deviation) including both experimental and simulation errors. Experimental data: (□) (49), (○) (50), (◊) (51), (▵) (52) and (+) (53). (Reproduced with permission. Ramos-Mendez J, Perl J, Schuemann J, McNamara A, Paganetti H, Faddegon B. Monte Carlo simulation of chemistry following radiolysis with TOPAS-nBio. Phys Med Biol. 2018; 63:105014.)

Independent reaction time. The pre-chemical stage from Geant4-DNA has been linked to the independent reaction time (IRT) approximation, which has been made available for access in TOPAS-nBio extensions via the TOPAS scoring manager. IRT (54, 55) is a computationally efficient alternative to the Brownian dynamics simulation for the calculation of chemical yields resulting from water radiolysis. IRT is a stochastic technique consisting of the sampling of reaction times and chemical reactions of pairs of species, independent from the surrounding neighbors. The method has been shown to be equivalent to full step-by-step Brownian dynamics (56). For a given pair of reactive species initially separated by a known distance, the aim is to determine at what time the pair reacts together. For that, the solution to the diffusion equation given by the Green function is inverted to retrieve the reaction time. The form of the Green function depends on the type of reaction (totally or partially diffusion controlled) and the charge between species. The position of reactive products was estimated using the position approach described by Clifford et al. (54) due to its ease in coding. The iterative sampling process is performed until a time cut specified by the user (e.g., the time to achieve the steady state of the system) is reached or all possible combinations between species are sampled.

Independent reaction time is computationally more efficient than propagating reactants using Brownian dynamics. However, as the number of species in the track increases, the efficiency decreases due to the iterative process. To mitigate this, Green et al. (55) proposed the sampling of reaction times only between chemical species separated by a distance of no more than d_max, the distance traveled by the most diffusive species before reacting (H₃O+ in this work). This is done in TOPAS-nBio. In addition, the particle track is binned spatially into a sparse three-dimensional matrix to further reduce computation time. For a given chemical species located at position , only those species contained in rectangular parallelepiped voxels fully or partially contained in a virtual sphere of radius d_max and centered at are used for the sampling. The voxel size is set to d_max/2 but it can be defined by the user.

The number of species and reactions currently implemented in the IRT of TOPAS-nBio are those used by Frongillo et al. (57) and recently by Plante and Devroye (56). This implementation includes a total of 72 reactions classified into six types covering first-order reactions, partially and totally diffusion controlled, between neutral and charged species.

Figure 6 shows G values as a function of time for 300 MeV (∼0.3 keV/µm) proton track segments of 100 µm in length calculated with both step-by-step Brownian dynamics using TOPAS-nBio with Geant4-DNA chemistry and with IRT from TOPAS-nBio. Statistical uncertainties in both calculations were below 0.5% (1 standard deviation). A revised physics list from Ramos-Méndez et al. (38) was used, where the relevant Geant4-DNA physics models included the elastic scattering from the CPA 100 implementation and the Bom model for ionization and excitation. Chemistry parameters with the electron thermalization distance were also from (38). As shown, the IRT and step-by-step results are in comparable agreement with the experiment. Both simulation methods used the same input data at the 1-ps starting point (within statistical uncertainties). Note that the different chemical species are tested for contact reactions at zero time by Geant4-DNA before IRT is performed. There is an indication that IRT has improved accuracy of yields at 1 µs. This is attributed to the more complete scheme of reactions and the use of partially diffusion-controlled reactions in the implemented IRT. The execution time was dramatically reduced with IRT, by a factor of ∼ 145, taking less than 6 s per proton track on a 2.7 GHz Intel^® Xeon^® processor.

FIG. 6.

G values as a function of time for Geant4-DNA/TOPAS-nBio step-by-step simulation (dash-dot line) and IRT (solid line). Experimental data from different sources are represented by symbols. For : (□) = ⁶⁰Co γ rays (58), (▪) = ∼2 MeV electrons (59), (•) = 20–22 MeV electrons (60), the latter with data scaled by a factor of 0.8 [see (61)] and plus sign (+) = 7 MeV electrons (62). For e^-_aq: (□) = ∼35 MeV electrons (63), (▪) = ∼45 MeV electrons (64), (○) = ∼40 MeV electrons (65, 66), (•) = ∼2.9 MeV electrons (67) and plus sign (+) = 20 MeV electrons (68). For H₃O⁺: (□) = 5 MeV electrons (69), (▪) = ⁶⁰Co and 8 MeV electrons (70), (○) = 3.5 MeV electrons (71) and (•) = 15 MeV electrons (72). For H₂O₂, “X” indicates ⁶⁰Co γ rays (58). For H₂: (▵) = ⁶⁰Co γ-rays (73). For : (▵) = ⁶⁰Co γ rays (74). [Reprinted (and expanded) with permission. Ramos-Mendez J, Perl J, Schuemann J, McNamara A, Paganetti H, Faddegon B. Monte Carlo simulation of chemistry following radiolysis with TOPAS-nBio. Phys Med Biol. 2018; 63:105014.]

Scorers

TOPAS provides flexibility in scoring quantities for further analysis. We extended the geometry-coupled scoring of TOPAS to cell and sub-cellular geometries. Note that the physical quantities one might want to score will depend on the scale and should be selected accordingly [e.g., macro-dosimetric LET vs. microdosimetric lineal energy (y) distributions].

Track-based scoring. While it is generally sufficient for macroscopic Monte Carlo simulations to provide volume and surface scorers, simulations at the sub-cellular scale often require scoring parameters based on each single particle track. For this purpose, the concept of n-tuples, (i.e., an ordered set of values, each containing n elements), has been incorporated. The contents of each different tuple were defined in a new extension class. We have provided a template for users to define the scored quantities, as well as a few options of n-tuple scoring. The n-tuples can be written out in files using ASCII, binary or root data formats. Two special cases of track-based scorers are described below.

SSB and DSB scorers. Simulations of radiation damage often focus on damage to the nucleus, in particular to the DNA double helix. Breaks in the DNA fall into two main categories: SSBs, i.e., breaks that involve only one strand of the DNA, and DSBs, i.e., breaks on both strands of the double helix within a certain defined distance, typically set to be within 10 basepairs. The complexity of the breaks, i.e., how many backbone sugars and bases are damaged within one damage site, is assumed to greatly influence the probability of repair. Accordingly, most studies define additional sub-categories such as a DSB+, a DSB with additional damage within a predefined surrounding basepair.

Defining an SSB and DSB scorer for generic representations of DNA is nearly impossible due to the inherent subtleties in the design of each DNA geometry. Such scorers need to know how basepairs are defined and how the geometries are set up within Geant4, e.g., as replicas, parameterizations or independent geometries. In TOPAS-nBio, these scorers are thus linked to a specific geometry. Currently, TOPAS-nBio offers two scorers that can generate geometry-dependent SSB and DSB distributions in a tuple format, TsScoreDSBFibre for the TsFibre geometry, and TsScoreDSBFractal for the TsFractal geometry. For both scorers, the user has to select an energy threshold to generate a break (e.g., 17.5 eV) and the number of basepairs between breaks that still count as a single break (e.g., 10 basepairs). These scorers can be used as templates to design scorers for other DNA geometries by updating the geometry dependence of the scorer for each specific implementation.

Another method to score DSBs without considering the actual geometry is one based on the density-based spatial clustering of applications with noise (DBSCAN) algorithm (75) that has been adopted by Geant4-DNA (76). DBSCAN is used to estimate the number of DSBs based on the distribution of energy deposition events in a given volume (e.g., nucleus) and a likelihood of clusters to form a DSB based on geometric considerations (e.g., the fraction of the volume covered by DNA). We have adopted the DBSCAN algorithm in TOPAS-nBio. A sample parameter file for this scorer is shown in Fig. 7.

FIG. 7.

Sample scoring parameter files. Left side: Scorer for the DBSCAN algorithm to determine DSBs without detailed geometry simulation. Right side: Parameter file to create a SDD scorer.

The standard for DNA damage scorer. TOPAS-nBio fully supports a new standard data format for DNA damage (SDD) scoring, which has recently been developed (39) both to score damage induction and as an input for repair modeling. The SDD scorer is a more generic version of a DNA damage scorer, such as the DSB/SSB scorers, and captures more information about radiation-induced changes to the DNA. However, similar to the SSB/DSB scorers, such a scorer requires detailed information about the simulated geometry. Accordingly, each TOPAS-nBio SDD-scorer has been linked to a specific DNA geometry class. We provide SDD scorers for the same two DNA geometries as for the DSB scorer (TsFibre and TsFractal). These scorers can again be used as templates to design scorer extensions for other DNA geometries. To record damage in an SDD format, users must select one of the two DNA structures, select the SDD fields to be filled and add the information for the DNA repair kinetic modeling that is required by the SDD format or that can optionally be filled. An example parameter section of an SDD scorer is shown in Fig. 7 for a nucleus containing the fractal walk DNA (Fig. 4B).

Graphical User Interface

The TOPAS text-based parameter system facilitates the setup of Monte Carlo simulations. To further ease the burden of setting up, running and reviewing results of TOPAS-nBio simulations such as virtual radiobiology experiments, we developed a GUI. The goal was to give users access to all the parameters available through the text-based parameter system while at the same time offering a visual interface to TOPAS-nBio (see Fig. 8). The GUI was implemented as part of TOPAS. Thus, while the GUI was designed for TOPAS-nBio users, it will eventually be included in standard TOPAS releases. While users will see the GUI as a whole new front end to TOPAS, the internal design is such that the GUI can be used to adjust values of TOPAS and TOPAS-nBio parameters and then run their simulation from those parameters, so that the already well tested TOPAS parameter control system is still in control. Users can also save the parameter setup with the GUI as a new set of parameter files to later run with or without the GUI. This allows users to run their revised simulations on non-GUI machines such as batch or cloud systems.

FIG. 8.

GUI developed for TOPAS-nBio. Changeable parameters can be set in the GUI, new geometrical components, sources, and scorers (shown) can be added via a drop box interface. The GUI can show the visualization of the simulated geometries and particle tracks once the “Run” button is pressed. Images can be exported to PDF.

The GUI allows users to set up complete simulations. A selection of predefined physics settings is adjusted to be optimal for the length scales of the simulation setup. The user may add geometric objects, such as cells, their nuclei, and their DNA, and scorers to the geometric components, such as a double-strand break (DSB) scorer or scoring DNA damage in the SDD format. In particular, all of the developments described here are available to view and change in the GUI. To allow such multiscale simulation design, the GUI provides a convenient method to zoom between the scales.

Comparisons and Validations of Track-Structure Simulations

Here we summarize our previously published studies validating the simulated track structure, comparing initial DNA damage induction to other codes and experimental data.

Stochastic proton interactions in fluorescent nuclear track detectors. In an effort to validate the simulated track structure, we compared the simulations to experiments performed with fluorescent nuclear track detectors (FNTDs) at the MGH proton beamline (77). Multiple measurements were taken at 2–4 mm intervals along pristine and spread-out Bragg peaks; however, due to limited microscope time, only a few data points were evaluated. We concluded that the intensity of the fluorescent signal observed with confocal microscopy should correspond to the number of secondary electrons stopping within a radius corresponding to the optical resolution of the microscope (207 nm). We compared the total intensity of the spot summed from the image data within our region of interest, defined by the microscope resolution and the simulated number ofelectrons stopping within the same radius of primary proton tracks. The resulting intensity plots at different depths were fitted using a skewed Gaussian. The comparison of the fit results of the Gaussian widths is shown in Fig. 9. We further found strong correlations between the FNTD track intensity and both the track-averaged LET and the frequency mean microdosimetric lineal energy y_F. While this study was limited by detector resolution, it provided a first validation of the Monte Carlo simulated track structure at the subcellular scale.

FIG. 9.

Validation studies for TOPAS-nBio. Panel A: Comparison of TOPAS-nBio track-structure simulations to measurements using FNTD detectors as described elsewhere (77). Shown are correlations of the width of skewed Gaussian fits to the simulated number of electrons stopping within a 207-nm radius of primary proton tracks and experimentally derived integrated brightness of FNTD tracks within a radius of 207 nm. Panel B: Comparison of the frequency of energy depositions in nucleosome-sized cylinders between different code's, TOPAS-nBio, MOCA8, OREC and CPA 100 used by Charlton et al. (78, 79). Panel C: Comparison of TOPAS-nBio SSB and DSB yields in dry plasmids and experimental data (79, 80). (Reprinted with permission. McNamara AL, Ramos-Méndez J, Perl J, Held K, Dominguez N, Moreno E, et al. Geometrical structures for radiation biology research as implemented in the TOPAS-nBio toolkit. Phys Med Biol. 2018; 63:175018.)

Simulations in a simple DNA geometry. To validate TOPAS-nBio using historical Monte Carlo track-structure data, a DNA segment using a representation developed by Charlton et al. (81) was modeled and simulated with both monoenergetic alpha particles and electrons to compare to the simulated results in (78). The energy deposited within the DNA segments was modeled using the Geant4-DNA physics processes and both the number and type of strand break within the segment were calculated and compared to that reported in Charlton et al. Breaks were defined when the energy deposited within the sub-volume (sugar-phosphate cylinder) was greater than 17.5 eV and then categorized into six different types: SSB, 2SSB, SSB+, DSB, DSB+ and DSB++. Differences between the Charlton et al. study and TOPAS-nBio included the generation of more low-energy events as well as more DSBs and fewer DSB++S. However, the overall occurrence of breaks showed similar trends to the Charlton et al. study [see (79) for details]. The observed differences were attributed to the physics models of the very low-energy events between the two track-structure codes, with TOPAS-nBio predicting more low-energy deposition events compared to the MOCA8 code, which was utilized in the Charlton et al. study. Better agreement was found between Geant4-DNA, CPA 100 and OREC simulations than with MOCA8, since the cross sections in MOCA8 are based on data for gaseous water, while the other codes are all based on calculations for liquid water (see Fig. 9B). This was consistent with that reported in another comparison study of dose point kernel distribution calculations (82).

Plasmid DNA. To validate TOPAS-nBio using experimental data we modeled the irradiation of a dry circular plasmid DNA model with 200 basepairs to compare to the measurements in Vysin et al. (80). The energy deposited within each component of the DNA plasmid (i.e., the sugar-phosphate backbone) was scored in an n-tuple scorer to determine the total number of SSBs and DSBs as a function of proton LET. Both the experimental and simulation data show an increase in both SSB and DSB yields with increasing LET. Both studies also show that protons are more effective at producing SSBs than DSBs in dry samples. The experimental data, however, predicted higher SSB to DSB ratios for all LET values by approximately a factor of 2 (see Fig. 9C). Differences between the two data sets are likely due to an oversimplified DNA model geometry or other experimental factors not accounted for in the simulation (79, 80).

Example of Complex Simulation with TOPAS-nBio

As a proof of concept, and to illustrate application of TOPAS-nBio to an artificial, complex simulation, we simulated protons of 3 and 50 MeV impinging on a chromatin section consisting of seven nucleosomes based on the TsFibre class, as developed in (26), to assess the difference between simulations with and without chemical reactions activated. Damage was scored on the backbones of the DNA geometry (see Fig. 10A). All physics interactions were scored as a damage to the backbone without requiring an energy threshold and a chemical damage was scored when a entered a backbone of the DNA geometry, resulting in the immediate removal of the from the simulation. Other radicals (e.g., solvated electrons) were tracked but not scored as damage inducing in this simulation unless they produced an . The propagation of chemical reactants may be even more important for larger volumes where DNA at more distant sites can be affected. Scavengers limit the distance of chemical propagation in the cell environment. The effect of scavengers was approximated in this demonstration by limiting the transport of chemical reactant propagation to a time of 1 ns. The propagation time can be adjusted through the parameter system to match the anticipated lifetime of chemical species in cell environments with different scavenger concentrations.

FIG. 10.

Simulation of a chromatid fiber irradiated with 3 and 50 MeV protons, including the chemistry of radiolysis. The center shows a representation of the simulated geometry and particle tracks, including the propagation of chemical species. Panels B–E and F–H: The results for 3 and 50 MeV proton irradiations, respectively. The plots demonstrate the effects of including chemical reactions. Here, only interactions with hydroxyl radicals are scored. Approximately one half of the DSBs include a reaction with hydroxyl, which not only influences the overall number but also the complexity of DNA damage (defined as number of lesions within a DSB).

Figure 10A shows the results of the simulations on the geometry displayed in that panel. Figure 10B shows the overall number of breaks (SSBs and DSBs) for 3 MeV proton irradiations and the number of breaks that involved interactions. Over one half of the DNA damage included at least one interaction, with some damage sites involving up to four interactions. Damage sites here are defined as a group of DNA lesions (physical energy depositions or interaction) with no more than 10 basepairs between neighboring lesions. Accordingly, damage sites can spread over significantly higher numbers of basepairs than 10; in our simulations, the largest DSB observed spread across 41 basepairs. Figure 10C (and F) shows the number of interactions involved in damage categorized as DSBs for 3 MeV (50 MeV) protons. It can be seen that the fraction of damages in DSBs is slightly less than for all breaks combined. DSBs are defined by having at least one damage pair on opposing strands within 10 basepairs.

Figure 10D and E (G and H) show the complexity of the DSBs, i.e., how many backbones were involved in the damage site for 3 MeV (50 MeV) protons. Multiple physical or chemical damage within the same backbone were counted as single instance. Figure 10D and G show damage including chemical reactions in dark brown and the complexity only from ionization events, i.e., damage sites not involving any chemical reactions in light beige color, and Fig. 10E and H show DSB complexity for damage sites caused entirely by chemical reactions. This latter damage has lower complexity (maximum 5) but would be absent in physics-only simulations. Furthermore, damage involving chemical reactions may be more difficult to repair, since attachments can alter the DNA molecules and structure and thus potentially the repair mechanism. Thus, while the damage patterns may be similar, repair probabilities between DSBs induced by physics events alone or a combination of physical and chemical events may differ. Multiple other factors such as chromatin structure (euchromatin vs. heterochromatin) also influence damage induction and repair (20, 83). In the presented simulations, only a small section of the TsFibre irradiated with a monoenergetic and monodirectional source was simulated. The propagation of chemical reactants may be even more important for larger volumes where DNA at more distant sites can be affected. The inclusion of scavengers or, as in this case, the limitation of chemical reactant propagation to a time of 1 ns, limits the distance of chemical propagation to more realistic cell environments. Overall, these simulations show that it is important to study differences in damage patterns and severity when chemistry is considered in order to better understand the processes within a cell.