Application of the Inverse Monte Carlo Method to High Energy Physics Data Analysis

1. Application of the Inverse Monte Carlo Methodto High Energy Physics Data Analysis William L. Dunn Quantum Research Services, Inc.* June, 2002 Effective 1 August 2002: Department of Mechanical and Nuclear Engineering Kansas State University

2. Preview • Nomenclature and Statement of the Problem • Inverse Monte Carlo Formalism • Prior Applications in Nuclear Engineering • Medical imaging • Radiative transfer • X-ray fluorescence • Wedge filter design • Potential Applications in HEP • Conclusions

3. General Monte Carlo Formulation Consider • where • y is a vector of independent variables • x is a vector of random variables • f is a joint probability density function (pdf) • X is the domain of x • K is a kernel or transformation function •  =  +  is a vector of m parameters

4. Monte Carlo as Quadrature Assume X, K, f , and  are known and exists for k = 1, 2,..., M. Form the quantity where the i are sampled from f*. Then, by the Strong Law of Large Numbers, and, by the Central Limit Theorem, we can estimate the uncertainty in the

5. The Inverse Problem • Direct Monte Carlo is just a form of quadrature, which is most advantageous for complex multidimensional problems • The standard application of Monte Carlo allows us to form estimates of expected values (definite integrals) • For inverse problems, where parameters in the integrand are unknown, the standard procedure is to assume values of these parameters, evaluate the integral(s), iterate until a pre-defined level of agreement between measured (or desired) values and Monte Carlo estimates is achieved; this is iterative in the simulation and thus inefficient

6. Inverse Monte Carlo Formulation If R has the known (measured or specified) values R* at M values of y and if  is unknown, then, for each yk; k = 1, 2,..., M, is a Fredholm Integral Equation of the First Kind, which defines a general class of inverse problems for the unknown 

7. Inverse Monte Carlo Solution Instead of iteratively applying direct Monte Carlo, formally apply Monte Carlo for a chosen pdf, f*, and score the results with the unknown a in the scores, viz. where the i are sampled from an arbitrary (but known) pdf, f*, defined and nonzero on X. Then forms a system of M algebraic equations in the m unknowns, .

8. Interpretation of IMC • The original integral equation with independent variable y has been converted into a system of algebraic equations at a specific set of yk which can be solved, in principle, by standard methods such as least-squares, maximum likelihood, successive approximation, or matrix inversion (if M=m and the dependence is linear). • In a sense, we have determined the “response” of the system for an arbitrary weighting of the unknown parameters, and we need only determine the proper weighting that best matches the solution to the known or measured responses.

9. Features of Inverse Monte Carlo • IMC is non-iterative in the simulation. • IMC applies, in principle, regardless of the dimen-sionality of the integral or the complexity of K or f. • The uncertainties in the recovered parameters can be estimated, using the Central Limit Theorem. • As in direct Monte Carlo, judgment, experience, and innovation in defining the model are rewarded. • On the other hand, it may be difficult to manipulate the unknowns within the Monte Carlo summation in order to form a manageable system of equations.

10. Applications in Nuclear Physics • Medical Imaging IMC model was tested with single-photon emission computed tomography (SPECT) phantom data to demonstrate scatter and attenuation compensation for a 2-D reconstruction. C.E. Floyd, R.J. Jaszczak, and R.E. Coleman, Inverse Monte Carlo: A Unified Reconstruction Algorithm for SPECT, IEEE Trans. Nucl. Sci.NS-32: 779 (1985). • Radiative Transfer Optical properties of inhomogeneous media can be recovered from surface measurements such as albedo and transmission. W.L. Dunn, Inverse Monte Carlo Solutions for Radiative Transfer in Inhomogeneous Media, J. Quant. Spectrosc. Radiat. Transfer29: 19 (1983). S. Subramaniam and M.P. Menguc, Solution of the Inverse Radiation Problem for Inhomogeneous and Anisotropically Scattering Media Using a Monte Carlo Technique, Int. J. Heat Mass Transfer34(1): 253 (1991).

11. Energy-Dispersive X-ray Fluorescence Ternary (Ni- Fe- Cr) sample 109Cd radioisotope ring source Shielded Si(Li) detector We seek to determine the weight fractions of all 3 elements in the sample from the total (primary and secondary) detected X-ray intensities.

12. Basic IMC Approach • Simulate to obtain the X-ray intensities incident on the detector for an assumed (reference) sample composition that contains all 3 elements. • Treat unknown composition as a perturbation about the reference composition. Thus, solve for where  is density w is weight fraction and subscript 0 refers to the reference sample. • Introduce the j in the scores by applying correction weight factors to force the interaction points in the unknown sample to be at the same locations as in the reference sample.

13. Intensity Models Using a second-order Taylor series expansion of the exponential terms, we obtain for k = 1, 2, 3, where A, B, C, D, E, F, G, and H are constants determined entirely by Monte Carlo, the Ik are measured total intensities, and the j are the unknowns to be determined.

14. Numerical Results No. Ni Fe Cr Max % Correct IMC Correct IMC Correct IMC error 1 10 10.1 80 80.2 10 9.7 3 2 20 20.3 60 59.5 20 19.8 -1.5 3 20 19.2 40 40.2 40 39.7 4 4 40 39.7 40 40.1 20 20.2 1 5 40 39.9 20 20.1 40 40.0 0.5 6 10 10.2 10 10.3 80 79.5 3 7 80 79.8 10 10.1 10 10.2 2 One Inverse Monte Carlo run for each sample using reference sample composition based on relative intensities. “Known” intensities estimated by direct Monte Carlo calculation.

15. Wedge Filter Design Photon beam modifiers are used in radiation therapy to flatten the intensity distribution and to create desired dose profiles within the patient. We consider the inverse problem of determining the wedge angles, , of lead filters placed between a 60Co source and a water phantom, given measured isodose curves within the phantom. In this case, the dose at position r in the phantom can be written

16. Geometry

17. IMC Model An IMC model can be written where and  are scatter and total attenuation coefficients, K is an unknown normalization factor, and Akji and Bkji are Monte Carlo scores for photons that do not interact in the filter and those that do, respectively. Rearranging, we obtain the matrix equation

18. Isodose Data For nominal 30° and 60° wedges,for which the isodose curves have relatively constant slopesof 30° and 60° over substantial distances in x (see notional Depiction at right), we took measured dose values at 60 positions, forming a measureddose matrix D for each wedge filter.

19. IMC Solution We input M=60 doses from the measured isodose curves, and solved the IMC model for wedge angle, obtaining Nominal Actual IMC Descriptor 30 7.2 7.1 60 16.1 16.9

20. Potential Application to High Energy Physics • Search for Particle Properties • Refinement of top quark mass • Higgs boson • SUSY models • Calorimeter Design • Electromagnetic calorimeter performance • Hadron calorimeter performance • Muon system design

21. Example: Top Quark Mass • Previously, Monte Carlo event samples generated for many discrete estimates of mt for a given channel (e.g., leptons+jets) and tested against data • The IMC approach uses single unknown of mt and generates events keeping track of their dependence around the assumed value; record dependence of mass on production cross section, top quark pT, decay rate into Wb, and kinematics of the particles produced. Then fit observed distribution (SET of jets) and fit to background plus signal. Invert to obtain mt that best fits model to the data. Single large simulation, with very good statistics

22. Conclusions • IMC provides a non-iterative simulation approach to solve a general class of inverse problems. • IMC solutions allow full phase-space simulation. • Each inverse problem is converted into a system of algebraic equations; the normal problems associated with existence and uniqueness of the solutions still apply. • IMC has been successfully applied to several problems in nuclear physics. • IMC offers promise for various problems in high energy physics. • The principal advantage is that simulation resources can be applied to obtain models whose coefficients have very good precision and that can be solved by algebraic techniques.