Statistics in AstronomyInitial question: How do you maximize the information you get from your data?Statistics is an art as well as a science, so it’s important to use it as you would any othertool: you have to know what you want to get out of your statistics before you can choosethe right way to apply them. Therefore, let’s start with some philosophy.First, statistics are not magic! If the data aren’t good enough, no amount of statisticalmanipulation will bring out the results. If your data are exceptional, the main results areoften evident even without lots of massaging. Nonetheless, there are plenty of cases in themiddle (where only good statistics will work), and it’s appropriate to be rigorous whenpossible anyway.Second, in astronomy systematics are often a major issue. An important principle isthat hypotheses are tested in bundles. When you make an observation, you aren’t measuringthe physical quantities of the system directly. You make implicit assumptions about yourknowledge of the instrument, and often need to make explicit assumptions about a basicmodel, to get the quantities that really interest you. Often there are enormous uncertaintiesas a result. This means that the statistical significance you quote may not be representative.For example, suppose I claim an astronomical result at 3 sigma significance. If I know myinstrument and ancillary assumptions perfectly, this will only happen by chance 3 times in1000, so that sounds good. But long and painful experience shows that systematics can eatthat up in a hurry, which is why people often adopt a much higher threshold for a significantresult. Example: when the third BATSE catalog of gamma-ray burst positions was released,the team had revised its position-finding algorithm. They said that it corresponded well tothe previous positions, because only 4% of the burst positions moved by more than 5σ. Now,5σ is a 3 ×10−7result, so obviously there were some significant systematics there!Third, when quoting significance you need to be careful about the number of trialsyou’ve performed, and in particular you need to be careful about hidden trials. Suppose I’mlooking for evidence of ESP by having subjects guess patterns on cards. At the end, I’mecstatic because there are more hits than expected, at the 95% significance level. Even if myexperiment is designed perfectly (and these usually aren’t), I have to realize that there maywell be 20 other labs doing the same thing, and only the positive results are quoted.So what’s a good approach? One thing you do not want to do is just use statistics blindly,any more than you would obliviously use the Saha equation for the center of our Sun (whereit doesn’t apply). An approach I find helpful is to try to think of how a certain statistical taskwould be done correctly (if only I had time and all the information I needed), then determinea reasonable approximation for my specific case. My upbringing as a postdoc, combined witha reasonable amount of experience, has made me comfortable with a generalized Bayesianapproach. Some of the characteristics of this approach are forward folding, the use of Poissonlikelihood, and a clear distinction between parameter estimation and model comparison.Let’s look at each of these in turn, then apply them.The idea of forward folding is that we should try to evaluate models based on directcomparison with the data. You may ask: isn’t this what we always do? Not necessarily. Sayyou want to estimate the temperature of a molecular cloud based on some observations. Howdo you do it? The way you might think of is, say, to get line ratios from your observations,then use an equilibrium model to get the temperature from the line ratios. This is backwardsfolding: take the data, then (figuratively) propagate it backwards to the molecular cloud.Forward folding would be having a model of the cloud including temperatures, densities, andso on, then a model of what emission this would produce, then a model of how that lightwould be received by your detector. You would then compare your “model data” with thereal data. I must say that in practice the two methods are often nearly the same, mainlybecause the detector is well-understood. However, for a given case you should at least givesome thought to whether some aspect of your detector is uncertain enough that backwardfolding could be a problem.The Poisson likelihood can be used any time your data come in discrete intervals (whichwe’ll call “counts”), and the counts are independent of each other. Schematically, we imaginedividing data space up into “bins”, which could be bins in energy channel of our detector,location on the sky, time of arrival, or any of a number of other things. Suppose that ina particular model m, you expect there to be micounts in bin i. Then if the model iscorrect the likelihood of actually observing dicounts in bin i of the data is, from the Poissondistribution,Li=mdiidi!e−mi. (1)Note that mican be any positive real number, whereas dimust be an integer. Note also thatthe sum of Lifrom di= 0 to ∞ is 1. The likelihood for the whole data set is the product ofthe likelihoods for each bin:L =Ymdiidi!e−mi. (2)This becomes better and better approximated by a Gaussian as miincreases.It’s interesting to talk to some people and ask why they did certain things with theirdata. For example, maybe they’ve taken a spectrum and put it into coarser bins. Whenasked why, they may respond “I put it into coarser bins to make the bin-by-bin statisticsbetter”. That is, if originally there were few counts per bin, they rebin so there are lots ofcounts per bin and thus they can use Gaussian statistics. Here’s a secret: when you makeyour bins coarser, you lose information, always! This has to be true, since by rebinning inthis way you no longer know where in the bin each count came from. If one uses Poissonlikelihoods, small numbers are fine. In fact, if you can manage, the best way to representyour data is to have bins so tiny that you expect either 0 or 1 count per bin. Then, you’regetting maximum information. There are times when your data don’t come in this count-by-count way. For example, when observing at Keck you don’t count individual photons. Inthat case, each bin of the data itself may have error bars associated with it. Then one canreplace the Poisson likelihood with a standard Gaussian; this is because, effectively, thereare many counts per bin.What