Adv. Space Res. Vol. 11, No.8, pp. (8)95—(8)99, 1991
0273—1177/91 $0.00 + 50 Copyright © 1991 COSPAR
Printed in Great Britain. All rights reserved.
KERNEL DENSITY ESTIMATORS APPLIED TO FAST TIMING HARD X-RAY OBSERVATIONS OF THE CRAB PULSAR I. R. Carstairs,* A. Bazzano,*** A. J. Court,* A. J. Dean,* N. A. Dipper,* M. J. Gorrod,* R. A. Lewis,* P. P. Maggioli,** F. Perotti,** E. Quadrini,** J. B. Stephent and P. Ubertini*** *physics Department, The University, Southampton, U.K **JFC, CNR, Milano, Italy ***JAS, CNR, Frascati, Italy tTESRE, CNR, Bologna, Italy
ABSTRACT The newly-proposed kernel density technique provides an objective way to search for pulsar signatures and to estimate the shape of the underlying source radiation emission. This new approach has been applied to some fast timing hard X-ray observations of the Crab pulsar obtained during the 1986 flight of the MIFRASO balloon-borne telescope, with a timing resolution of 0.33msec over the energy range 15 to 300 keV. The resulting pulse profile shows the classical bimodal density seen in other observations, and is compared with the traditional binned histogram method. The associated H,~test for pulsar periodicity searches is reviewed as a flexible alternative to conventional tests. INTRODUCTION The pulsar in the Crab Nebula is one of the best studied and most familiar objects in modern astronomy, yet many questions remain unanswered. A particular emphasis is needed for fast timing hard X-ray observations of the Crab, as this region is closely related to the most energetic phenomena occurring at the central pulsar. Studies of the pulse shape can reveal much information on the geometry, characteristics and mechanisms of the emission regions /1,2/. THE TELESCOPE A detailed description of the MIFRASO (MIlano, FRAscati, SOuthampton) hard X-ray telescope (15-300 keY) may be found elsewhere /3/. The main photon detector system is designed around a modular construction of 6mm thick NaI(Tl) scintillation counters (area 2700 sq.cm) and a high pressure xenon gas proportional counter (area 900 sq.cm). The proportional counters provide higher spectral resolution (typically 12% FWHM at 60 keV) for brighter sources than the scintillation counters (25% at 60 keV). The time of arrival of incident photons is recorded with a precision referenced to universal time of typically 0.5 ms (dependent on the telemetry bit rate). Directionality is achieved by means of a passive collimator system providing a hexagonal aperture of 2.6 degrees FWHM. The X-ray telescope is mounted in an appropriate structure which includes an elevation drive system and azimuthal control. The microprocessor based steering system was designed to provide a positioning accuracy and stability better than 0.1 degrees. Observations are possible in either drift scan or tracking modes. DATA ANALYSIS The arrival time of each accepted photon event was transformed to the corresponding barycentric arrival time using the JPL DE200/LE200 ephemerides /4/ with measurements of the Crab pulsar period provided by contemporary radio observations /5/. The pulse profile was then determined by folding the events at a single period interpolated from the radio ephemeris.
(8)95
(8)96
I. R. Carstairs etal.
13000
•
12000
1~1
12000
I
r
000000 955000 950000 ~045000 ~ 040000
E.FI
•[i~
12100
J~L~ 910000
090000
12000
885000~~ 0
0.2
0.4
0.9 0.8 1.0 1.2 relôtlvi phasi
1.4
1.8
1.0
2.0
0.4
0:6 ~ ~.0 L.2 r.I~Llva ph~..
~
.0
Figure 1: The pulse profile of one 1986 Crab nebula observation as obtained using the KDE and binned histogram techniques Various techniques have been proposed for estimating unknown density distributions. Traditionally the histogram estimator has been used to bin the event phases where all information within a bin is lost. A new method has been employed to extract the pulse profile of the source; the Kernel Density Estimator (KDE) technique /6/. This relies on the assumption that any periodic function can be estimated by using a kernel (or weighted) function evaluated at each phase point, together with a harmonic series of Rayleigh moments of the phases. The kernel function employed was that suggested by Swanepoel /7/ and the number of harmonics was found by applying the method suggested by Hart /8/. This new method allows an objective estimate of the shape of the pulsed emission to be made, from which to derive a more realistic ~description of the underlying source emission mechanisms. Figure 1 illustrates the pulse profile found using KDE and a binned histogram. The profiles show similarities with previous hard X-ray measurements: two narrow pulses separated by a phase of approximately 0.4, with a pronounced interpulse structure between the two. The secondary pulse also shows some sign of a broad shoulder on the leading edge as seen previously, and the relative heights of the two peaks are also consistent with earlier results /9/. STATISTICAL TESTS Tests for pulsed emission in X and 7-ray astronomy fall into two groups: parametric tests, where tile shape of the periodic signal is assumed, and non-parametric tests, where nothing about the shape of the signal is assumed. A periodicity search using the latter class of test is conducted by folding the event arrival times with an appropriate trial period and possible period derivative to form a set of phases. The resultant phase density distribution, an estimation of which forms the ‘light curve’ or ‘pulse profile’, is then tested against the null hypothesis of phase uniformity, 7-h,. If ‘i~is rejected with a confidence greater than some chosen significance level, the alternative hypothesis, fli, is accepted: that there is a periodic signal present within the data. Care must be exercised in distinguishing a true source from an artifact of the experiment, for example, deadtime or bitrate effects. Comparison with an identical periodicity search for an equal number of off-source event times is helpful in distinguishing tile latter effect. There is a wide variety of statistical procedures with which to test 7I~against The problem with most tests is that they are powerful in detecting only certain kinds of phase density distributions, whilst other distributions may pass unnoticed. The two statistical tests most commonly used are the Pearson and Rayleigh tests for uniformity: ~
Kernal Density Estimators 12
12
—,—~ i
~~~
I
~l
T
I
~
HioH :
13)’
—
•,
(8)97
/
IPS
7~
6-
.
~
l.r)nl period
I~~I~~II_
I
.~°
LrIol pr~rincl / It’S
Figure 2: Periodicity search using a) Rayleigh and H,1, tests, b) Pearson test with 20 bins and phase origin o = 0.0. THE PEARSON TEST Tile Pearson test /10/ defines an origin for the phase domain and divides it into a number of equal bins. The test is a (smoothed) non-parametric test for uniformity because no assumption is made concerning the underlying form of the phase density distribution with the number of bins, b, acting as a smoothing parameter. The test is easy to invoke but retains an intrinsic uncertainty in its significance through arbitrary choice of bin size and phase origin as a narrow peak may be split between bins, reducing the power of the test. The statistic is not invariant under a global phase change of origin o, and narrow features may be ‘massaged’ into one bin, subjectively increasing the power of the test. Finally the test is powerful in detecting ~ frequency oscillations which, for large b, are unlikely to correspond to a physical source, thus 7~t~ may be falsely rejected. The Pearson test is powerful in detecting narrow pulses of duty cycle .~. and is thus popular for the analysis of Crab or Vela-type light curves with 20-40 bins per phase domain. THE RAYLEIGH TEST The Rayleigh statistic /11/ is a parametric test for power in the fundamental harmonic of the trial period, although it may be used non-parametrically, for example, in Fast Fourier Transform techniques. If tile phase density distribution contains narrow peaks the sensitivity of the test is lowered as significant power is contained in higher harmonics away from the fundamental. In such cases other tests such as Pearson’s are more powerful. There is one important case where the Rayleigh test fails: if the density distribution is bimodal, as seen in the double-peaked Crab and Vela pulsar light curves, moments cancel out as the phase angle between 0,and the the trigonometric test may falsely reject Ni. the peaks approaches ~r The Rayleigh test was first advocated for the use in YHE 7-ray astronomy /12/ as ‘the uniformly most powerful test for the uniformity of a circular (phase) distribution’. This implies that the test has the largest power of all possible tests for all possible types of phase density distribution. As noted above, there are circumstances under which the Rayleigh test may falsely reject fl~and because one test is UMP it does not follow that another is not. It may be shown through the Neyman-Pearson lemma that both the Pearson and Rayleigh tests are both UMP under certain conditions /13,14/. Tile Rayleigh test is powerful in detecting broad, unimodal phase density distributions such as those seen in the high energy emission of accreting X-ray binary systems. THE H, 1, TEST The bias in most conventional tests towards particular phase density distributions, and the sub-
(8)98
i. R. Carstairs cc aL
jective choice in applying smoothed tests, have motivated the search for a powerful, flexible and objective statistical test for uniformity. Such a test should be nearly as powerful as Pearson’s for detecting narrow peaks and nearly as powerful as the Rayleigh test for broad peaks. This required flexibility suggests a smoothed test, but the smoothing should be chosen in an objective and optimal manner from the data itself. A test which fulfills these requirements has recently been proposed /15/ and is called the Hth test. The Hth test extends the Z~,,test /16/, summing the first in harmonics up to a truncation point in = ñz. For m = 1 the test reduces to the Rayleigh test, powerful for detecting broad peaks, whilst a choice of m ~ 10 makes the test powerful in detecting narrow peaks. However, a subjective choice of m remains unsatisfactory and is decided using a method suggested by Hart /8/. The use of th is a maximum likelihood estimation of m /17/ where the likelihood function is easy to compute from the data itself and provides an optimal choice for in with firm statistical justification. In the analysis of X-ray observations with large signal-to-noise ratios and sample sizes the Hart function is fairly flat and poor at discriminating the choice of ñl, however, the choice of ñz will remain objective and optimal. With m = th the test statistic changes because ñ~iis now a random variable. The behaviour of Hth was investigated through simulations for 1 < m < 20 and the probability distribution parameterized /15/. This test is also complimentary to the kernel estimation of the pulsar light curve. Figure 2. show a comparison between the above tests on a subset of the data. Numbers in parentheses refer to the optimal number of harmonics selected by the Hart rule. The vertical dashed line corresponds to the interpolated contemporary radio period and the error bar for the Pearson test indicates variations in detection significance as the phase origin is slid across the bin width for that trial period. The powerful and objective behaviour of the H,t, test confirmed using these observations suggests that it is a good test for pulsar periodicity searches. The value determined by the tests also confirms the value used for folding the event times to within uncertainties associated with the timing and barycentring routines. It should be noted that folding after a period search, to investigate the pulse profile from a source with unknown period, the number of independent trials, together with possible oversampling within the independent period space, will degrade the statistical significance. CONCLUSIONS The performance of the newly-proposed Hm test for uniformity has been compared to two traditional tests using a periodicity search on a small set of hard X-ray observations of the Crab pulsar. This new test performs as powerfully as the most appropriate traditional tests, with no subjective choices or variations in significance level. Its use should be encouraged for a wide variety of different sources and instruments in order to investigate its behaviour more fully. REFERENCES 1. 2. 3.
4. 5. 6. 7. 8. 9. 10. 11.
K.S. Cheng, C. Ho, M. Ruderman, Astrophys. J. 300, 500 (1986). K.S. Cheng, C. Ho, M. Ruderman, Astrophys. J. 300, 522 (1986). R.E. Baker, G. Barbaglia, A. Bazzano, L. Boccaccini, A. Bussini, A. Carzaniga, A. Court, A.J. Dean, N.A. Dipper, G. Ferrandi, H. Haskell, C. Lapadula, R.A. Lewis, D. Maccagni, M. Mastropietro, R. Patriarca, F. Perotti, V.F. Polcaro, E. Quadrini, D. Ramsden, S. Sembay, R. Spicer, P. Ubertini, G. Villa, D. Whatley, Nuci. Instr. Meths. 228, 183 (1984). E.M. Standish, Astron. Astrophy. 114, 297 (1982). A.G. Lyne, R.S. Pritchard, .Iodrell Bank Crab Pulsar Timing Results, Monthly ephemeris published by University of Manchester Nuffield Radio Astronomical Laboratories. O.C. de Jager, J.W. Swanepoel, B.C. Raubenheimer, Astron. Astrophys. 170, 187 (1986). J.W.H. Swanepoel, Commun. Stat. 16(6), 1835 (1987). J.D. Hart, J. Stat. Comp. 21, 95 (1985). R.B. Wilson, G.J. Fishman, Astrophy8. J. 269, 273 (1983). K. Pearson, Nature, 72, 294 (1905). J.W.S. Rayleigh, Theory of Sound, McMillan and Co. 1894, p.35.
Kernal Density Estimators
12.
13. 14. 15. 16.
17.
(8)99
A.I. Gibson, A.B. Harrison, I.W. Kirkman, A.P. Lotts, 3.H. Macrae, K.J. Orford, K.E. Turver, M. Walmsley in: Proc. ml. Workshop on V.H.E. 7 Ray Astronomy, ed. P.V. Ramana Murthy, T.C. Weekes, Ootacamund, India, Smithsonian Inst., U.S.A., p.97 (1982). G.A. Korn, T.M. Korn, Mathematical Handbook for Scientists and Engineers, New York: McGraw-Hill (1968). P.G. Hoel, Introduction to Mathematical Statistics, New York: Wiley (1971). O.C. de Jager, B.C. Raubenheimer, A.R. North, H.I. Nel, G van Urk, Astrophys. J. 329, 831 (1988). R. Buccheri, K. Bennett, G.F. Bignami, J.B.G.M. Bloemen, V. Boriakoff, P.A. Caraveo, W. Hermsen, G. Kanbach, R.N. Manchester, J.L. Masnou, H.A. Mayer-Hasselwander, M.E. Ozel, J.A. Paul, B. Sacco, L. Scarsi, A.W. Strong, Astron. Astrophys., 128, 245 (1983). J.H. Friedman, in: Data Analysis in Astronomy II, ed. V. di Gesü,L. Scarsi, P. Crane, .
J.H. Friedman, S. Levialdi, New York and London: Plenum Press, p.3-15 (1986).