Linear Algebra Applications in Neutrino Experiments
Linear Algebra plays a crucial role in various neutrino experiments, enabling solutions in weighted least squares, energy estimation in EXO-200, and signal processing in LArTPC detectors. From Cartesian coordinates to minimizing uncertainties, linear algebra techniques contribute significantly to data analysis in neutrino physics research.
Download Presentation
Please find below an Image/Link to download the presentation.
The content on the website is provided AS IS for your information and personal use only. It may not be sold, licensed, or shared on other websites without obtaining consent from the author. Download presentation by click this link. If you encounter any issues during the download, it is possible that the publisher has removed the file from their server.
E N D
Presentation Transcript
Review of Linear Algebra Applications in Some Neutrino Experiments Xin Qian BNL 1
Introduction Linear Algebra (LA) has a very long history: First appears in The Nine Chapters on the Mathematical Art Systematically introduced by Rene Descartes An ancient Chinese book Application of LA is very broad: We selected a few examples in the field of neutrino physics to illustrate its power Cartesian coordinates in geometry 2
Weighted Least Squares (WLS) ( ) ( ) T = R S R S 1 T M S R M C M Minimizing T leads to a solution with Linear Algebra of : (vector) measurement : ( : response matrix connecting signal to measurement : Covariance matrix describing uncertainties ( ) 1 unknown vector ) signal = T S C R C M S and ( ) 1 = 1 T C R C R S C A. C. Aitken Proc. R. Soc. Edinburgh 55, 42 (1935) An approximation or a special case of the maximum likelihood approach With unbiased measurements, this solution is the best linear unbiased estimator (BLUE) 3
Improved Energy Estimator in EXO-200 Search for 0 with 136Xe During the initial operation, observed correlated noise in the APD light readout Liquid Xe Time Projection Chamber Light waveform of an example noise event An optimal energy estimator based on linear algebra to reduce correlated noise: JINST 11, P07015 (2016) 4 Nature 510, 229 (2014)
Improved Energy Estimator in EXO-200 ( ) = R E R E + 1 ( ) ( ) T T ( ) T N M N M N N C N T = R E R E 1 T M R E C M C M : ( unknown nuisance vecto ) excess noise r : (vector) APDs' waveform : (vector) response : ( : Covariance matrix measure d wit minimizing T = = ed by the biconjugate gradient stabilized alg. (BCGSTAB: a numerical method) R E + + T T T D = 1 D N 0 R C M M R C R S R + N u n no n c k w s al a ) energy r ( ) 5 ~ 1.43 1 I C Eqs. s l o v h no is e e v ent s Major challenge: large dimension: 2048 (time bins) x 70 (APDs) ~ 1.43x105 difficult to calculate (RTC-1R)-1 5 More details in JINST 11, P07015
Signal Processing in Single-Phase LArTPC Collection Y unipolar shape Induction (U,V) bipolar shape Number of ionized electrons Field Response Signal on Wire Plane Electronic Response Signal to be digitized by ADC (Signal Processing) Number of ionized electrons time Goal of signal processing is to reconstruct the distribution of ionization electrons from TPC waveforms JINST 13 P07006/7 (2018) 6
Matrix Inversion and Deconvolution ( ) ( ) Tensor R is Toeplitz Response only depends on the relative positions between signal and measurement T = R S R S T M M M : (matrix in wire and time) TPC waveform : ( electron distribution : (tensor ) r e s p R unk n o w n m t i a r ) ionization x S (RTR)-1 can be achieved through the deconvolution technique with fast Fourier transformation (FFT) O(N2.4~3) O(N logN) One of top 10 algs. in 20th century o ns e Immense dimension of M and S: 2400 (chs.) x 8256 (time bins) ~ 2 x 107 Extreme challenge in matrix inversion! 1D deconvolution SP in LArTPC described in B. Baller JINST 12, P07010 7
Deconvolution Technique Time domain = ) ( ) S t dt M( ) ( M(t) t R t t S(t) 0 0 t Fourier transformation M = R( ) S( ) ( ) Frequency domain M( ) S( ) = ( ) F R( ) PMT Signal Processing in Daya Bay: NIMA, 895, 48 Back to time domain Inverse Fourier transformation S(t) Wiener filter: optimized by power of signal and noise in freq. domain N. Wiener (1964) 8
TPC Signal Processing: 2D Deconvolution ( ) 0 0 0 1 0 1 M ( ) ( ) ( ) R (t t) ( ) ... i i i t M R S + = + + = + + dt t R t t S t S t + 2-D 1-D Measured waveform deconvolution ( ) ( ) ( ) R ( ) S ( ) ... deconvolution 0 1 1 i i i R1 represents the induced signal from (i+1)th wire signal to ith wire ( ) ( ) ( ) ( ) ... ( ) ( ) ( ) ( ) ... ( ) ( ) ... ... ... ... ... ( ) ( ) ( ) ( ) ( ) ( ) ... ( n S M M R R R R R R R R S S 1 0 1 2 1 1 n n 2 ... 1 0 3 ... ( ) ( ) 2 ... ( ) ( ) 2 n n = ( ) ) M M R R R R R R R R S 1 ( ) 2 3 0 1 1 ( ) n n n n 1 2 1 0 n n n The inversion of matrix R can again be done with deconvolution 2-D Fast Fourier Transformation JINST 13 P07006/7 (2018) 9
Intermediate Summary So far, problems with N N measurement unknown i) Direct matrix inversion, ii) Numerical solution, iii) Fast Fourier Transformation, What about N << N ? measurement u nknown ) ( ) ( 1 1 = 1 T T Original solution fails, as S R C R R C M 1 T the matrix can not be inverted! R C R 10
Solution: Compressed Sensing A powerful mathematical technique to recover sparse signal from under- determined (linear) system ||S||0 counts the number of non-zero entries in S = R S Sovle M L0: minimize ||S|| , subject to M n n S M = R S Find the most sparse solution (NP-hard) 0 Breakthrough: L0 (combinatorial) solution can be approximated by the (iterative linear) L1 regularization (E. Candes, J. Romberg, T. Tao arXiv- math/0503066) minimize T ( ) ( ) L0 reg. ! L1 reg. O N O m N ( ) 2 = R S + | | M S i i So far, one of the top 10 algs. in this century 2log ( = + or generalize to T , ) | | L M S S i i Minimizing through coordinate descent (an iterative approach), J. Friedman et al. J. stat. soft. 33, 1 Allow positivity constraints. Implemented as RESS package 11
Wire-Cell Tomographic Event Reconstruction Three-dimensional Imaging for Large LArTPCs , XQ, C. Zhang, B. Viren, M. Diwan, JINST 13, P05032 (2018) Sparsity, positivity, and connectivity information are added through compressed sensing (L1 regularization) ~ large charge true hits ~ zero charge fake hits 12
Wire-Cell on the MicroBooNE Data 2015, after several hours of CPU running Bee link (upper) 1D deconvolution + L0 regularization 2018, 1 min of CPU running Bee link (lower) 2D deconvolution + L1 regularization JINST 13, P05032 (2018) 13
Neutrino Detection in MicroBooNE Steps Steps Steps Neutrino: Cosmic Neutrino: Cosmic Neutrino: Cosmic Rejection Power Rejection Power Rejection Power Selecting neutrinos (or rejecting cosmic muon rays) on a surface- running LArTPC is a big challenge! Hardware trigger Hardware trigger Hardware trigger 1:20k (~26 per event) 1:20k (~26 per event) 1:20k (~26 per event) 1 1 1 Software trigger Software trigger 1:800 1:800 0.04 0.04 Offline light filter 1:200 0.01 14
Compressed Sensing in Selecting Neutrinos 4.8 ms drift window Coin. with beam 20-30 TPC activities 40-50 PMT activities LArTPC is a slow detector, running on surface leads to significant challenges in selecting neutrinos The key is to associate the correct T0 (light info.) to the corresponding TPC cluster Combinatorial problem compressed sensing! 15
4.8 ms drift window Matching Principle 20-30 TPC activities All possible hypotheses One flash many or zero TPC clusters within corresponding active volume (activities in inactive volume) One cluster at most one flash (inefficiency in the light system) Light signal proportional to (reconstructed 3D) charge Known light acceptance given position Predicted vs. Measured light pattern with Compressed Sensing 16
Core Charge-Light Matching Algorithm Overall test statistics to be minimized = + + + 2 2 ij 2 p 2 p 2 p : Measured Light Pattern : Predicted Light Pattern : Uncertainty : ith Light Flash : jth PMT : kth Charge Cluster k M P 1 2 3 , i j 2 Comparison of the measured and predicted light pattern b M M a P ij ik ikj ij ij = k 2 ij 2 ij M i 2 j Each Charge Cluster can only be used once 1 a ik = k 2 p 1 2 0.01 b i 2 ij Observed light flash may not correspond to any charge cluster = 2 p Aggressively pursue charge- light matching Additional cuts to examine the light mismatch events 2 2 0.025 , i j = 2 p | | a Compressed sensing to select the best pairs 3 ik , i k https://lar.bnl.gov/wire-cell/ 17
LArTPC Charge-Light Matching After Charge-Light Matching Top view Side view Front view After 3D Imaging Steps Neutrino : Cosmic Muons Offline light filter (late light correction) 1 : 200 Many-to-many charge-light Matching 1 : 7 18
Through-going-muon rejection and effective detector boundary Vertical Direction (Y) Anode E field (X) Cathode measurement Vertical Direction (Y) prediction Beam Direction (Z) Vertical Direction (Y) Drift direction (X) After correcting t0 for all cosmic muons in one event Steps Neutrino : Cosmic Muons Many-to-many charge-light matching 1 : 7 After through-going-muon rejection 1:1 19
What about N Compressed Sensing What else? Iterative Approach << N ? measurement unknown 20
Iterative LA Approaches A famous example is Bayesian unfolding G. Agostini, NIMA, 487, 362 (1995) & arXiv:1010.0632 Observation Truth ( ) ( ) n i + R S 1 th n + ji = ( i 1) n iteration S M j j ( ) n l R R S li jl l l : (vector) measurement : ( : smearing matrix M S R Reconstructed unknown vector ) signal to be unfolded N N signal measurement Positivity, but no sparsity Used in many experiments extracting -A cross sections Equivalent to the more general Maximum Likelihood Expectation Maximization (ML-EM): A. P. Dempster et al. J. Royal Stat. Soc. B39, 1 (1977) L. A. Shepp and Y. Vardi, IEEE Trans. Med. Imaging 2, 113 (1982) NEXT 0 : JINST 12, P08009 21
Iterative Approaches with Linear Algebra Iterative approach based on Linear Algebra to approximate nonlinear fits E.g. Track trajectory and dQ/dx fit in MicroBooNE BCGSTAB & Coordinate Descent Iterated Weighted Least-Squares fit to reduce bias in the Neyman or Pearson 2 in counting experiments arXiv:1807.07911 and earlier works M 2 Poisson = + 2 ln M M ( ) 2 M 2 Neyman = M ( ) 2 ( ) : solved in ith iteration : to be solved n + N 2 M + n = 1 i 2 i 2 Pearson = + 1 1 i 22
Combined Neyman Pearson (CNP) Chisquare ( ) 2 M 2 3 2 Poisson 2 Neyman A linear combination of Neyman and Pearson 2, CNP 2 leads to a better approximation of Poisson 2, advantages in computation in certain scenarios 2 M M ( ) 2 1 3 2 Poisson 2 Pearson 2 M 1 3 ( ) X. Ji, W. Gu, XQ, H. Wei, C. Zhang arXiv:1903.07185 = 2 Neyman + 2 Pearson 2 CNP 2 23
Linear Algebra in Designing Experiments M. Wilkings et al. arXiv:1412.3086 Conditional Covariance Matrix, also Gaussian Process arXiv:1709.05681 24
Summary Linear Algebra is a very powerful tool in experiments Least square approach: : N N constraint i) Matrix Inversion (modest N) ii) Numerical Solution (large N) iii) Fast Fourier Trans. (Toeplitz) N N unknown Applications Selected Examples Signal Processing MicroBooNE TPC and Daya Bay PMT waveform analysis Data unfolding EXO-200, Wiener-SVD Event Reconstruction Wire-Cell 3D imaging, Charge-light matching in LArTPC : constraint i) Compressed Sensing L1 reg. (Sparsity and Posit ii) ML-EM iterative approach (Prior and Positivity info.) unknown Iterative approach: ivity info.) Applications Selected Examples Numerical solutions Coordinate Descent, BCGSTAB Data unfolding ML-EM (Bayesian unfold) in NEXT, Xs Reducing bias Iterative Least Weighted Squares Approx. NL fit Trajectory & dQ/dx fit in LArTPC Also experiment s designs 25
Stopped Muon Rejection Drifting Drifting Beam Vertical Vertical Search for the Bragg peak in dQ/dx near the stopping point is the key, also important for PID Determination of track direction is crucial to reject stopped muons 27
Trajectory and dQ/dx fit in Single-phase LArTPC U view Measurement: Three 2D projections (wire and time) with charge information Broadening caused by diffusion, detector response, signal processing V view Hypothesis is a 3D thin line (many points) with charges (3+1) x N unknowns for N 3D points A non-linear fit approximate with two linear fits: trajectory and dQ/dx Allow further iterations to refine the solution W view 28
( ) = + + + Overall Test Statistics , , , T x y z Q T T T T j j j j U V W reg dQ/dx fit: Trajectory fit: Unknowns Measurements 2 2 i q ( ) 2 = T / / dis U V W q R Q / / U V W i Uij j 2 i ij q j j i = , T U 2 i q : pixel in 2D projection : 3D trajectory point i j i U T = , : smearing coefficients calculated with a fixed trajectory R ( ) ( ) ( ) ( ) 2 2 Uij = + 2 i 2 2 ( ) , , , , dis U U U U x y z x t t x y z j i j j j j i j j j j : bin size in U view, :bin siz in drift time t U x 29
Stopped Muon Rejection Drifting Drifting Bragg peak Beam Vertical Vertical Steps Neutrino : Cosmic Muons After through-going-muon rejection 1:1 London Wenqiang Gu Cooper-Troendle After stopped muon rejection 7.6:1 30
Generic Neutrino Detection Efficiency ~ 90% for charge-current (CC) interaction Selected CC interaction: ~ 11.5 k with 5e19 POT vs. 4.3 k in PRL 123, 131801 (2019) @ 50% purity High-performance generic neutrino detection achieved 31 Hanyu Wei
Remaining Challenges in Signal Processing JINST 13 P07006/7 (2018) Still worse performance for the induction wire plane for prolonged topology (long in drift time) because of the signal ROI finding 32 Hanyu Wei Brooke Russell
ROI finding as Image Segmentation UNet: arXiv:1505.04597 Time Tick Deconvolved with filters: field, and electronic effects removed Channels Time Tick ROI finding - Segmentation with UNet 33 Channels Haiwang Yu
Deep Learning-Signal Processing Performance on ProtoDUNE Data Current Signal Processing Deep-Learning Signal Processing Time Ticks Channels Significant Improvements in prolonged tracks 34 Haiwang Yu
Wiener-SVD Unfolding Approach LArTPC Signal Processing SVD Unfolding Approach Statistical (and systematic) uncertainties and smearing of finite response function leads to unstable unfolded results Unfold the ionization electron distribution from measured TPC waveforms High-frequency electronic noise and the finite resolution of detector response leads to unstable deconvolved results Regularizations (e.g. constraining smoothness) is added to stabilize results (limiting variance) Regularization Bias Wiener filter is added to stabilize the results, minimizing Optimize the strength of regularization to balance between bias and variance Strictly speaking, should also optimize the choice of binning (multi-dimensional scan) ( ) 2 ( ) ( ) ( ) ( ) E F M R S exp exp: Expected signal : Filter F S 35
Wiener-SVD Unfolding Approach When the expected signal is known , Wiener filter can be used in SVD unfolding to obtain an optimized result in one shot Choose binning < smearing resolution Apply Wiener filter : smearing matrix : Wiener filter : smearing in the effective frequency domain : expected signal in effective freq. domain i s n effective freq. domain 2 Ui 1 2 M N + ( ) s ( ) = log log L s 2 i 2 i 2 i d + s : signal to be unfolded : Likeliood function L s = U D V = s T R W i 2 i 2 i 2 i n d s ( ) R D W 2 : Expected "noise" in the effective frequency domain : Expected signal in the effective frequency domain N M 2 i 1: (white) noise i n Naturally balances bias (sys.) and variance (stat.) a smaller total uncertainty with no scanning of regularization strength More details in X. Li et al. JINST 12, P10002 36 Wei Tang Xiaoyue Li
An Example of Energy Spectrum Unfolding JINST 12, P10002 ( ) : regularization strength : bias 2 = + MSE bias Variance b Wiener SVD unfold outperforms the standard SVD unfold approach with a simplified procedure 37 Wei Tang Xiaoyue Li
Linear Algebra in Experiment Design 100% cancellation of flux uncertainty with one reactor, one near and one far detector Prog. Part. Nucl. Phys. 83, 1 Daya Bay ~95% 38