KAHIBARO
41.6. Combining Results
Table of Contents
ROOT output
On high performance computing systems you will often run many GATE jobs in parallel, each producing one or more ROOT files. Before combining them, confirm that every job used the same configuration except for elements that are meant to vary, such as random seeds or beam angles in a parameter study. The ROOT structure, for example tree names and branch names, must be identical in all files you plan to merge.
For purely additive quantities such as hits, singles, coincidences, or event-level phase space, combining ROOT output is usually a matter of concatenating trees. In ROOT this is often done with the hadd utility, which merges trees with identical structure. When using hadd, you provide an output ROOT filename and a list of input files, and the resulting file contains all entries from the input trees, one after another. This is appropriate if each file represents a disjoint subset of events from the same conceptual simulation. It is not appropriate if the trees differ in structure or if the same events were partially written to several files.
When combining with Python tools, for example with uproot, you read the same tree and branches from each file and then concatenate them in memory before saving to a new ROOT file or to another format. Make sure that any event identifiers or run identifiers are treated carefully. For many analyses it is enough to keep them distinct per file, but if you rely on unique event IDs across the full dataset you may need to remap IDs, for example by adding an offset depending on the job index. The important point is that the combined file should be logically equivalent to one very long single-run simulation.
Histograms that have already been filled by actors or analysis codes can also be combined by summing bin contents. When summing histograms from several runs, you should also combine their statistical uncertainties. If $H_i$ is the content of bin $i$ in a given run, and $\sigma_i$ is the bin uncertainty, then for a sum over runs indexed by $k$ you have a total bin content $H_i^{\text{tot}} = \sum_k H_{i,k}$ and a total uncertainty
For independent runs, combine per-bin uncertainties in quadrature:
$$
\sigma_i^{\text{tot}} = \sqrt{\sum_k \sigma_{i,k}^2}.
$$
This rule holds for count-based quantities where bins from different runs are statistically independent. If a bin contains simple counts, the uncertainty for each run is typically $\sigma_{i,k} = \sqrt{H_{i,k}}$ under Poisson statistics.
Sometimes you want to average several simulations instead of summing them, for example to estimate a mean response over several independent seeds. In that case you compute the mean per bin and its uncertainty using standard statistical formulas for the mean and standard error, not by simply dividing the summed histogram by the number of runs without adjusting the error bars. For any combined ROOT result, record which files were merged and how the merge was performed, so that you can reproduce or modify the combination at a later stage.
Dose images
Dose actors often produce voxelized dose outputs in image formats such as MHD or NIfTI. When running on an HPC system, a common pattern is to split the total planned number of primary events into multiple jobs, each writing its own dose image. Since each job simulates an independent subset of histories on the same geometry and scoring grid, you can combine the images voxel by voxel.
Before combining dose images, verify that all jobs used the same grid definition. This means the same image dimensions, voxel spacing, origin, and orientation. If any of these differ, you cannot safely add images voxelwise without resampling, which belongs to image processing and not basic result combination. Also confirm that all images represent the same physical quantity, such as dose in Gray per run, or total dose in Gray for each partial exposure.
Most dose actors can optionally output both dose and dose squared, and often also track the number of hits or histories per voxel. These additional outputs are essential if you want to compute combined statistical uncertainties. The basic combination rule for independent simulations is that the dose in each voxel is additive. If $D_{v,k}$ is the dose in voxel $v$ from run $k$, then the total dose is
For independent partial histories on the same grid, total dose per voxel is:
$$
D_v^{\text{tot}} = \sum_k D_{v,k}.
$$
If your runs all simulated the same number of primary events, and you wish to estimate the average dose per run, you can also compute the voxelwise mean dose by dividing $D_v^{\text{tot}}$ by the number of runs. This is mainly useful for convergence studies. For normal clinical or research applications, users typically care about the total dose from the full set of primaries, so they keep the summed image.
When dose squared is available, you can compute voxel uncertainties across multiple runs. Let $D_{v,k}$ and $D^2_{v,k}$ be the dose and dose squared for run $k$ in voxel $v$. The total dose and total dose squared are the sums over runs. The variance of dose in voxel $v$ can then be estimated from the combined quantities, using standard relations between the mean of a quantity and the mean of its square. From this variance, you obtain the standard deviation and relative uncertainty for each voxel. This provides a more accurate picture of the Monte Carlo noise than simply examining a single run.
If your workflow uses images normalized by total simulated activity or total delivered monitor units, you must be careful during combination. Either ensure that each partial run uses the same normalization, in which case voxelwise summation remains valid, or accumulate unnormalized quantities and apply the final normalization only after combining. Once you have a combined dose image, you can feed it into downstream tools for visualization or dosimetric analysis, treating it as if it came from a single large run.
Statistics
Combining results from multiple HPC jobs is not only about summing data, but also about correctly propagating statistical information. For Monte Carlo simulations, each independent job provides an estimate of some quantity with an associated statistical uncertainty. When you combine $N$ independent jobs, you improve precision, and the overall uncertainty decreases approximately as $1/\sqrt{N}$ if all jobs simulate the same number of events.
For count-based quantities, such as the total number of detected photons or coincidences, Poisson statistics are a good approximation. If job $k$ yields a count $C_k$, then the total count is $C_{\text{tot}} = \sum_k C_k$ and the combined standard deviation is
For independent Poisson counts:
$$
\sigma_{\text{tot}} = \sqrt{\sum_k C_k} = \sqrt{C_{\text{tot}}}.
$$
From this you can compute the relative uncertainty as $\sigma_{\text{tot}} / C_{\text{tot}}$. This pattern extends naturally to histogram bin contents, where each bin behaves like a count.
For mean quantities, such as the mean energy deposited per event or the mean depth of interaction, each job produces a sample mean and a sample variance based on its subset of events. To combine them, you use formulas for the pooled mean and pooled variance, which take into account both the per-job means and the number of events in each job. If each job has the same number of events, the combined mean is simply the average of the job means. The combined standard error then follows from the pooled variance divided by the total number of events. If jobs differ in event count, weight each job by its number of events when computing the overall mean and variance.
Dose statistics require special attention because dose is often estimated per voxel. If you only have a single dose value per voxel from each job, and no dose squared, you can still use the fact that independent runs approximate independent Gaussian estimates for large numbers of events. For a given voxel, its combined value is the sum of the per-job doses and the uncertainty approximately scales as $1/\sqrt{N}$ if all runs are identical in size. However, this only provides a rough scaling, not a precise voxelwise variance. If you anticipate combining many jobs, it is better to configure your dose actors to produce enough information, for example dose and dose squared, so that you can compute a rigorous combined uncertainty map.
When planning large HPC campaigns, decide in advance which statistical quantities you will need. If you want to report confidence intervals, convergence tests, or comparisons to experimental data, design your actor outputs accordingly. After combination, always check simple global statistics, such as total number of events, mean dose in a uniform phantom, or integrated counts, to ensure that the combined result matches expectations given the sum of events from all jobs. This final check helps you detect missing files, double counting, or misaligned normalizations before you move to detailed scientific analysis.
Views: 12