{'fraction_non_alphanumeric': 0.0315471618149551, 'fraction_numerical': 0.012528805657828516, 'mean_word_length': 4.6210079654582, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Indirect Dark Matter Detection from Dwarf Satellites: Joint Expectations from Astrophysics and Supersymmetry

We present a general methodology for determining the gamma-ray flux from annihilation of dark matter particles in Milky Way satellite galaxies, focusing on two promising satellites as examples: Segue 1 and Draco. We use the SuperBayeS code to explore the best-fitting regions of the Constrained Minimal Supersymmetric Standard Model (CMSSM) parameter space, and an independent MCMC analysis of the dark matter halo properties of the satellites using published radial velocities. We present a formalism for determining the boost from halo substructure in these galaxies and show that its value depends strongly on the extrapolation of the concentration-mass (c(M)) relation for CDM subhalos down to the minimum possible mass. We show that the preferred region for this minimum halo mass within the CMSSM with neutralino dark matter is ~10^-9-10^-6 solar masses. For the boost model where the observed power-law c(M) relation is extrapolated down to the minimum halo mass we find average boosts of about 20, while the Bullock et al (2001) c(M) model results in boosts of order unity. We estimate that for the power-law c(M) boost model and photon energies greater than a GeV, the Fermi space-telescope has about 20% chance of detecting a dark matter annihilation signal from Draco with signal-to-noise greater than 3 after about 5 years of observation.


Abstract.
We present a general methodology for determining the gamma-ray flux from annihilation of dark matter particles in Milky Way satellite galaxies, focusing on two promising satellites as examples: Segue 1 and Draco. We use the SuperBayeS code to explore the best-fitting regions of the Constrained Minimal Supersymmetric Standard Model (CMSSM) parameter space, and an independent MCMC analysis of the dark matter halo properties of the satellites using published radial velocities. We present a formalism for determining the boost from halo substructure in these galaxies and show that its value depends strongly on the extrapolation of the concentration-mass (c(M)) relation for CDM subhalos down to the minimum possible mass. We show that the preferred region for this minimum halo mass within the CMSSM with neutralino dark matter is ∼ 10 −9 − 10 −6 M ⊙ . For the boost model where the observed powerlaw c(M) relation is extrapolated down to the minimum halo mass we find average boosts of about 20, while the Bullock et al (2001) c(M) model results in boosts of order unity. We estimate that for the power-law c(M) boost model and photon energies greater than a GeV, the Fermi space-telescope has about 20% chance of detecting a dark matter annihilation signal from Draco with signal-to-noise greater than 3 after about 5 years of observation.

Introduction
If the dark matter in the Universe is comprised of stable, weakly-interacting massive particles (WIMPs), in many instances this leads to the prediction that WIMPs will self-annihilate into Standard Model (SM) particles that may be visible with the upcoming generation of high-energy particle detectors. If high energy gamma-rays are produced, there are several promising sources within our own Galactic environment where the annihilation radiation from WIMPs may be visible. A detection of annihilation products from multiple sources, in possible concert with detections from colliders and underground labs, will be required to conclusively establish the nature of the dark matter in the Universe.
Each source of annihilation radiation has its advantages and disadvantages. Because of its close proximity and high dark matter density, the flux from annihilation radiation will be the largest in the direction of the Galactic center. However, uncertainties in the empirical determination of the central density profile [1,2,3,4] and in the contamination from gammaray sources that are not of dark matter origin [5] may hinder the extraction of a dark matter signal from this region. Both of these systematics may be somewhat alleviated by searching for annihilation radiation at a few degrees offset from the Galactic center [6], but even in this region a full analysis and understanding of the spectrum of the astrophysical backgrounds is required.
The high mass-to-light ratios and the relative proximity of the dwarf spheroidal satellite galaxies (dSphs) of the Milky Way make them also excellent independent targets that have been widely considered in the literature [7,8,9,10,11,12,13,14,15,16]. Their status as potential targets for indirect detection has become even more interesting recently, given that the known number of satellites has more than doubled in the past few years [17,18,19,20,21], coupled with the discovery that all of these satellites share a common dark matter mass scale within their central 300 pc [22]. While the actual signal of gammaray flux from the Milky Way dSphs is smaller than the flux from the Galactic center, the astrophysical gamma-ray backgrounds tend to be reduced in the direction of these objects, as most of them are located at high Galactic latitudes. Additionally, these dSphs have low intrinsic emission from astrophysical gamma-ray sources; not only do they have mass-to-light ratios greater than ∼ 100 in most cases, but all of them within ∼ 400 kpc have strong upper limits on HI gas content [23].
The astrophysical contribution to the calculation of the annihilation flux from any source can generally be divided into two components: one component arising from the smooth halo contribution, which is proportional to the density squared distribution in the halo, and one component arising from the flux due to bound substructures within each of these halos. While it has long been recognized that the presence of substructure in dark matter halos can have a significant effect on the annihilation rate of dark matter particles [24], theoretical calculations of this boost factor from substructure have varied by orders of magnitude because of large uncertainties in both the density profile of substructures and their distribution within the parent halos [7, 10, 25? ]. Numerical simulations are now reaching the necessary resolution to characterize the density distribution of substructures [26,27], resolving substructure down to four levels in the hierarchy [27]. These, the highest resolution simulations, model the evolution of the dark matter and do not include baryons. In principle, simulations of this kind can be used to provide accurate estimates for the expected boost factor only in objects with negligible baryonic mass fractions (e.g., dwarf galaxies). However, in regions like the Galactic center where baryons play an important role in the dynamics, baryonic effects, such as encounters of substructures with stars in the Galactic disk or bulge [28], or the backreaction of the dark matter distribution in response to disk formation [29], must be included. Future simulations that include more physics and even higher resolution will be required to fully characterize boost factors in general circumstances.
It is clear that, in order to extract an unambiguous detection of dark matter annihilation radiation, a full understanding of all astrophysical uncertainties is required. With a goal of characterizing these uncertainties, in this paper we present an algorithm for the calculation of the gamma-ray flux from dwarf satellites, accounting for both uncertainties in the smooth halo distribution and the halo substructure distribution. We introduce a method for scanning the parameter space and determining the best fitting dark matter distributions from the kinematics of stars in these satellites using a Monte Carlo Markov Chain (MCMC) analysis. We combine the parameter constraints from the satellite stellar kinematics with the constraints on the parameters of the Constrained Minimal Supersymmetric Standard Model (CMSSM), with a goal of obtaining predictions accounting in a realistic way for all relevant sources of uncertainty of the flux from the dwarf satellites. The predicted regions that we delineate will provide guidance for future gamma-rays experiments for testing the predictions of neutralino dark matter in the CMSSM self-consistently within the context of CDM. As examples, we apply our algorithm to two particular satellites, both of which are known to be strongly dark matter dominated: the classical satellite Draco, which is located at 80 kpc from the Sun, and Segue 1, which is a newly-discovered satellite at 23 kpc from the Sun. We provide these two as examples and leave the work of ranking all the satellites in terms of their expected gamma-ray luminosity to a future paper.
As a part of our analysis, we provide a new analytic solution to the equation governing the boost factor from halo substructure. Our solution is particularly useful because it allows for mass functions and halo concentrations to be free functions of host halo mass that can be implemented at each level in the substructure hierarchy. This is particularly important when implementing results from recent numerical simulations which show that the normalization of the mass function is reduced by up to 50% at the next level of hierarchy [27]. Using this solution, we show that the uncertainty in the boost factor is dominated by the extrapolation of the dark matter halo concentration versus mass relation down to mass scales that are currently unable to be resolved in CDM simulations. Assuming an optimistic power law extrapolation, we find mean boost factors ∼ 20, in agreement with recent numerical extrapolations for Milky Way mass halos [27]. Assuming a concentration-mass relation that is linked to the small-scale power spectrum as in the model Bullock et al. [30], however, leads to boost factors of order unity. As an additional application of this analytic formula, we solve for the minimum dark matter halo mass at each point in the CMSSM parameter space, and find that the typical range for the minimum mass CDM halo is 10 −9 − 10 −6 M ⊙ . This result updates and generalizes the previous calculations of the minimum mass halo in the context of Supersymmetric CDM [31,32,33,34,35].
The paper is organized as follows. In Section 2 we review the formalism for determining the gamma-ray flux from dark matter annihilations. In Section 3, we review our assumptions for the CMSSM and the SuperBayeS software that scans the CMSSM parameter space. In Section 4 we review the formalism to determine the best fitting dwarf satellite dark matter halo profiles, under the assumption of the CDM model. In Section 5 we present our calculation of the probability distribution for the boost factor and the resulting differences relative to the smooth halo flux. In Section 6 we discuss some detection prospects for present observatories, and in Section 7 we present our conclusions. Our main conclusions for the detectability of the flux are summarized in Figure 11.

Annihilation Flux and Gamma-ray spectrum
Following standard methods [1,36] the gamma-ray flux from particle annihilations can be derived via where ℓ is the line-of-sight distance, ℓ ± = d cos θ ± r 2 t,DM − d 2 sin 2 θ, d is the distance to the galaxy, θ is the line-of-sight angle from the center of the galaxy, < σν > is the average cross section for annihilation, r t,DM is the tidal radius of the dark matter halo, and m χ is the WIMP mass. Here, is the number of photons above energy E produced per annihilation and the resolution window function is For Fermi, θ res is approximately 10 arcminutes at the energies we consider. It is clear from Equation (1) that an accurate prediction for the flux entails incorporating relevant uncertainties both for the astrophysical quantities and for the particle physics model. In the present analysis we constrain the density of the dark matter halo, ρ(ℓ), from the kinematics of the stars, while the annihilation cross section, WIMP mass and annihilation channels we use are derived in the following section.

The Constrained Minimal Supersymmetric Standard Model
For the calculations in this paper, we will assume that the dark matter particle consists of the lightest stable supersymmetric particle. More specifically, we focus on the case of the CMSSM. In this section, we review the basic definitions and parameters of the CMSSM, discussing specifically the SuperBayeS code that we use to explore the CMSSM parameter space and how this feeds into our calculation of the flux.

Neutralino Dark Matter in the CMSSM
Supersymmetry (SUSY) provides a compelling and well-motivated extension to the Standard Model that naturally contains a dark matter candidate [36]. Supersymmetry postulates a symmetry between bosons and fermions -every boson (e.g. gauge bosons) has a fermonic partner (e.g. gauginos) and every fermion (leptons, quarks) has a bosonic partner (sleptons, squarks). The particularly well-studied R-parity conserving, weak-scale softly broken SUSY models provide both natural solutions to the "fine tuning problem" and a natural dark matter candidate [37]. The former is achieved through cancellation of quadratic divergences of one loop quantum corrections to the Higgs mass. The latter is due to a conserved discrete symmetry (R-parity) that prohibits the lightest SUSY particle from decaying to only SM particles.
Enlarging the particle sector in this manner greatly increases the number of free parameters that specify the model; even the most minimal form of SUSY (MSSM) introduces over a hundred new parameters. Such a large number of free parameters makes the efficient exploration of the MSSM parameter space very challenging. The naive method of exploring the likelihood surface on a regularly-spaced grid is clearly inadequate, as the required computational effort scales exponentially with the number of dimensions considered. Furthermore, present-day constraints on SUSY phenomenology are fairly indirect and do not allow for meaningful constraints on models with so many degrees of freedom. A popular and well motivated simplification is achieved by demanding SUSY parameter unification at GUT (Grand Unified Theory) scales. This so-called constrained MSSM (CMSSM) now limits the number of SUSY parameters to 4 continuous and 1 discrete parameter [38]: the common gaugino mass, m 1/2 , common mass for scalars, m 0 , trilinear scalar couplings, A 0 , (all of which are specified at the GUT scale, M GUT ≃ 2 × 10 16 GeV) the ratio of the higgs expection values, tan β, and the sign of the "µ term", sgn(µ) [37,38]. We shall denote the CMSSM parameters as It has been recently demonstrated [39,40,41,42,43] that the value of some SM parameters can strongly affect the predictions for some of the observable quantities, in particular the relic neutralino abundance (see Fig. 4 in Roszkowski et al. [41]). Therefore, it is not sufficient to specify the values of C and fix the relevant SM parameters to their current best-fit values, but the latter must be introduced as "nuisance parameters" and marginalized over, in order to account for their impact in the predictions. The most relevant SM parameters are namely, the top quark mass, the bottom quark mass, the strong coupling constant and electromagnetic coupling constant, respectively. In the context of the CMSSM, specification of the parameter set (C , N ) allows for derivation of a full suite of predictions for low-energy observables.
The package SuperBayeS , developed by Ruiz de Austri et al. [40] and Trotta et al. [43] embeds several codes in a MCMC framework to derive from (C , N ) the SUSY mass spectrum (via SoftSusy [44]), the neutralino relic abundance (via DarkSusy [45] or MicrOMEGAs [46,47]), SUSY corrections to various Higgs sector quantities (employing FeynHiggs [48,49,50,51]) and branching ratios of rare processes (using Bdecay [52]). The CMSSM and SM parameter space is explored by SuperBayeS using an MCMC Metropolis-Hastings algorithm, or, more recently, by employing the more efficient and robust "nested sampling" algorithm [43,53,54]. The parameters are then constrained by applying all available constraints on the low-energy observables, including the WMAP 5-years determination of the relic abundance, sparticle and Higgs mass limits, branching ratios of rare processes, electroweak observables and direct constraints on the SM quantities (for a detailed discussion of the likelihood, see [40,41,43]).

Priors in the CMSSM
The final outcome of the CMSSM analysis is a list of samples drawn from the posterior distribution, P(C , N |D) obtained via Bayes' theorem: where D denotes the combined data described above, L(C , N ) ≡ P(D|C , N ) is the likelihood function and P(C , N ) is the prior distribution for the CMSSM and SM parameters. From the posterior one can then derive the probability distribution for any function of the quantities (C , N ) one is interested in, for example the neutralino-proton interaction crosssection (relevant for direct dark matter detection experiments, see [55]), the gamma-ray and antimatter flux from the galactic center (relevant for indirect detection searches [56]) and Higgs-sector physics (of interest for Higgs-boson searches [57]). The posterior in Equation (6) should be dominated by the likelihood, L, so that the prior influences vanish for strongly constraining data (for more details, see e.g. [58]). However, it has been found that this is not currently the case for the CMSSM -namely, the available data are not sufficiently constraining to determine the CMSSM parameters in a prior-independent way (see Trotta et al. [43] for a detailed analysis). This means that some of the posterior constraints on C are somewhat dependent on the chosen prior distribution. The fundamental reason for this is that the mapping between high-energy CMSSM parameters and low energy observable quantities is highly non-linear, due to the nature of the Renormalization Group Equations, and therefore even fairly strong low-energy constraints are relatively mildly informative about the quantities one is interested in, namely C . It is however expected that this issue will be resolved once the LHC will deliver direct observations of the SUSY mass spectrum [59].

CMSSM samples and derived quantities
For the results in this paper, we use a nested sampling chain for the parameter space spanned by (C , N ), containing approximately 45, 000 samples. We assume throughout a positive sgn(µ) (motivated by consistency with the measured anomalous magnetic moment of the muon). We adopt the chain resulting from the analysis in Trotta et al. [43], with a flat prior with the limits denoted in table 1 (the "4 TeV" limits in Ruiz de Austri et al. [40]).
From those CMSSM chains, we derive for each sample the value of the WIMP mass, m χ , its annihilation cross section, < σν >, and number of photons produced in the annihilation above 1 GeV, N γ (1 GeV). As discussed below, we choose 1 GeV because it provides a conservative lower bound for the expected signal energy window from CMSSM dark matter annihilation. Furthermore, we compute the value of the minimum mass halo, m min , as explained in section 5.4 below.

Dwarf Satellite Kinematics
We follow the standard formalism for analyzing stellar line-of-sight velocities from dwarf satellites. In this section, we review the relevant formulae so as to establish notation and conventions that we use throughout this paper.

Theoretical Modeling
We consider the satellites as two-component systems consisting of stars and dark matter (e.g. [60,61]). The potential is assumed to be spherically symmetric, and the system is taken to be completely pressure-supported (no rotational support). This is seen to be a good description of many of the dwarf satellites [62,63,64]. With these assumptions the Jeans equation is where ρ ⋆ is the stellar density, σ r is the stellar radial velocity dispersion, β ≡ 1 − σ 2 t /σ 2 r is the velocity anisotropy parameter, and σ t is the stellar tangential velocity dispersion. The mass of the system, M(r), is defined as the total dynamical mass, which is the sum of the contributions from the stars and the dark matter. The line-of-sight velocity dispersion is where R is the projected radius onto the sky and I ⋆ (R) is the stellar surface density. The use of subscript "th" will become apparent below when Equation (8) is fed into our statistical analysis. Observationally, the integrated masses within ∼ 300 pc for all of the Milky Way satellites are very similar, consistent with ∼ 10 7 M ⊙ independent of the dwarf galaxy luminosity [22]. To first order, this fact simplifies the selection of the best targets for flux detection to those that are closest to the Sun [65]. However, as we discuss below, including the effects of a prior relation for the maximum circular velocity distribution for dark matter halos somewhat modifies this simple estimate. In general, the primary factors that determine the best sources are those objects that have the best signal-to-noise ratio, accounting for the astrophysical backgrounds. Some often-discussed targets include Segue 1 (23 kpc), Ursa Major II (32 kpc), Willman 1 (38 kpc), and Coma Berenices (44 kpc) (e.g. [65,66,67]) all of which were discovered since 2004. The half-light radii of these objects are ∼ 10-100 pc, and given their high velocity dispersions of ∼ 4 − 6 km/s [68], these objects are consistent with being darkmatter dominated dSphs [22,67,68]. The nearest of the more well-known (classical) dSphs include Ursa Minor (66 kpc) and Draco (80 kpc); previous calculations show that these two objects have similar predicted fluxes [10].
Since our main goal in this paper is to discuss the methodology for robustly predicting fluxes and including the boost in the calculation of the gamma-ray flux, we restrict our analysis to two example dSphs: Draco and Segue 1. The kinematic constraints we derive for Draco are much stronger, as this object has been well-studied both from the standpoint of its photometry and kinematics. Segue 1 is a more-recently discovered satellite that appears to be strongly dominated by dark matter and the least luminous galaxy known [67]. However, there are only 24 stars with measured velocities in Segue 1 [67], and, as we show below, the errors on its mass and flux are much larger than the respective values for Draco. The surface densities of these objects are fit by King [69] and Plummer [70] profiles, respectively. For Draco, the King core radius is r king = 0.18 kpc, and the King tidal radius is r t = 0.93 kpc [71]. For Segue 1, the one-component Plummer fit gives a Plummer core radius of r pl = 35 pc [72]. In the Plummer profile, ρ ⋆ falls off as 1/r 5 in the outer regions, so there is no natural definition of the stellar tidal radius, in contrast with the King profile. In this case, we conservatively assume that the stellar tidal radius is given by the position of the outermost observed star, which is located at a projected radius of R = 50 pc [67]. In principle, the stellar surface density parameters could also be estimated jointly with our other model parameters, however in the present work we choose to fix them instead to their best fit values as given above.
To model the density profile of the respective dark matter halos of these objects, we use an Einasto profile, which is defined as The Einasto profile has been shown to be a good fit to CDM halos with the Einasto index, n, ranging from ∼ 3 − 7 [73,74,75]. For our purposes, this profile is also convenient for two separate reasons: (1) the profile has a well-defined mass, which is important when we calculate substructure boost factors below, and (2) the profile does not diverge towards the center of the halo, which is convenient when calculating the gamma-ray flux. As presented in Equation (9), there are three parameters in the Einasto profile that we must determine from the data: the log-slope index n, the scale radius r s , and the scale density, ρ s . It will be convenient in the following to replace the Einasto density and radius variables (ρ s , r s ) with the implied halo maximum circular velocity and radius of maximum circular velocity (V max , r max ). The dark matter halo density profile is then specified in terms of the parameter set The final quantity we must specify is the velocity anisotropy, which enters both directly in the Jeans equation in Equation (7) and in the equation that relates the observed line-ofsight velocity dispersion to the underlying stellar radial velocity dispersion Equation (8). This quantity is unconstrained by line-of-sight velocity data [76], so in order to allow for general models, we model the velocity anisotropy as with four free parameters, We note that this parametrization is slightly more general than that used in Strigari et al. [61], allowing here for the power law index η in addition to the anisotropy scale radius, r β , and the asymptotic inner and outer anisotropies, β 0 and β ∞ . It should also be noted that there is an intrinsic degeneracy between the logarithmic slope of the density profile in Equation (9) and anisotropy. Even in the simplified case of constant anisotropy, this degeneracy restricts how well the central slope n of the halo may be contrained from stellar kinematics [76].

Likelihood function of dSph model parameters
Above we have discussed our theoretical modeling; we now turn to our description of the line-of-sight velocity data [67,77]. We begin by noting that the line-of-sight velocities from dSphs are well-described by Gaussian distributions [62]. The observed velocity distribution is a convolution of the intrinsic velocity distribution, arising from the distribution function, and the measurement uncertainty from an individual star. The probability of obtaining a set of data d given a set of model parameters H , V is described by the likelihood [10] where H are the parameters that describe the density profile of dark matter and V are the stellar velocity anisotropy parameters. The product is over the set of n stars, and u is the bulk velocity of the galaxy in the direction of the observer. As expected, the total error at a projected position is a sum in quadrature of the theoretical intrinsic dispersion, σ th,i (H , V ), and the measurement error σ m,i . Often, kinematic data from dSphs is presented in terms of the velocity dispersion in bins of projected radius. In these cases, it is useful to have an expression for the likelihood function similar to Equation (13) that is free of any terms associated with the measured velocities of individual stars. An expression of this form can be found by replacing u and σ by their respective maximum likelihood values,û andσ, where the latter quantities are obtained from a standard maximum likelihood procedure using Equation (13). Proceeding with this approximation, and also neglecting the measurement uncertainty in comparison with the intrinsic dispersion (which is a good approximation for the bright satellites), the likelihood function in Equation (13) can be reduced to , and ln(V max ) (blue, dotted). The prior distributions are truncated at V max = 3km s −1 as described in the text. Increasing negative powers of V max causes the posterior to be more "biased" toward lower mass solutions. As a result, the posterior corresponding to these different priors differ. This is now a product over the number of bins in projected radius, N, each with a velocity dispersionσ b . The number of stars in a given bin is N b . (Note the difference in the normalization relative to the expression presented in Strigari et al. [65], due to a typographical error in Strigari et al. [65]). Note that the quantity σ th,i = σ th,i (H , V ) is computed via Equation (8). As we did above in the context of the CMSSM parameters, we now go over to the posterior probability distribution function (pdf) for the parameters of interest, which is again given by Bayes' theorem where P(H ), P(V ) are the prior pdf's for the halo and velocity anisotropy parameters, respectively, which we take here to be uncorrelated. We will deal with the issue of priors in more detail below. The task is now to explore numerically the posterior, Equation (15), in order to determine credible regions on the model parameters. In previous similar work involving parameter estimation from dwarf satellites [61,65], the likelihood function was directly integrated over the range of model parameter space. This method was accurate but could be time-consuming, particularly in the case of large numbers of parameters to marginalize over. Rather than direct numerical integrations, in this paper we explore the posterior distributions using Markov Chain Monte Carlo techniques. Before discussing our MCMC methodology, we now turn to the discussion of the priors entering Equation (15).

Priors for dSph parameters
Choosing a prior in accordance with the physical situation and our degree of prior beliefs on it is an important aspect of Bayesian analysis, as variations in the priors can lead to sizable , and ln(V max ) (blue, dotted). The prior distributions are truncated at V max = 3km s −1 as described in the text. Note that the trend in mass within 300 pc with prior reverses here compared to the case in Figure 1. This traces back in part to the fact that the mass is best constrained near twice the half-light radius [78]. For Draco, 300 pc is within the half-light radius and when the CDM r max -V max prior is imposed, lower V max halos are forced, on average, to be more concentrated at 300 pc. In Figure 1, 300 pc is beyond the half-light radius of Segue 1, and priors that favor larger V max give larger extrapolated masses. differences in the posterior whenever the likelihood is not very strongly peaked. We account for this prior information in our notation by the inclusion of the appropriate infinitesimals. For example, a "uniform" prior probability in x will be denoted as P(x) = d(x) whereas a uniform prior probability in ln(x) is represented as Regardless of these definitions, if data strongly constrains some parameter, or a given combination of parameters, then the prior information should not have much bearing on the result, as the posterior is dominated by the likelihood. However as we will see for Segue 1, with only 24 radial velocities spanning only ∼ 50 pc in radius, priors can have a significant effect.
Guidance as to how to choose the priors for our dSph model parameters can be gleaned from cold dark matter simulations which give precise predictions of halo abundances given a halo shape and mass. We note however, that these simulations do not provide the probability of observing a galaxy, which itself depends not only on dark matter physics, but also on star formation and baryonic physics. But if we assume that primarily the halo mass, and not shape, affects stellar physics we may draw the halo profile prior probabilities directly from simulations, with some additional simple inputs to account for gas physics.
As discussed above, we describe the halo in terms of the parameters n, r max , V max . We model the r max prior (conditional on the value of V max ) as a log-normal distribution, which provides a good description of the r max − V max relation measured in the Via Lactea [79] and Aquarius [27] simulations over the entire range of V max . From the Aquarius simulations [27], we estimate this relation, given a V max , to be (16) where σ log(r max ) = 0.22 is a conservative scatter in the Aquarius subhalos for the entire range of V max . There are no statistics published (as of this paper) for the parameter n in the Einasto profile, though from Springel et al. [27] it is reasonable to assume a uniform prior in 1/n [i.e.
The choice of prior for V max is a more delicate issue, for a couple of reasons. One issue relates to the probability that a given CDM halo has a particular value of V max , and a second issue relates to how well the line-of-sight velocity data itself constrain particular values of V max . On the former point, V max is the primary parameter that relates to the astrophysical process of low-mass galaxy formation: small galaxies with more shallow potential wells are expected to have low star formation rates, so the actual V max prior is expected to be more shallow than that inferred from the CDM subhalo mass function, max of Springel et al. [27]). On the latter point, from the perspective of the line-of-sight data, large V max values are not constrained by the data. Given the form of the r max prior, large V max values correspond to large r max values that may fall outside the extent of the stellar profile. Therefore, these values cannot be directly constrained by data and accordingly become dominated by the prior. Figure 1 and Figure 2 exemplifies this behavior. The left panel of Figure 1 shows the posterior constraints on the mass of Segue 1 within 30 pc with four different prior assumptions: , and d(ln V max ). The right panel shows the mass within 300 pc using the same respective priors. Prior choice has little effect on M(30 pc) since this value is well constrained by line-of-sight velocity data. In contrast the M(300 pc) posterior is dominated by prior choice, as the radius of 300 pc lies outside the measured stellar distribution. In contrast Draco, whose stellar extent extends beyond 300 pc, M(300 pc) is well constrained. We note that the prior behavior in Figure 1 and Figure 2 does not contradict the results presented in Strigari et al. [22], where a uniform prior in ln[M(300pc)] was taken for the entire satellite population, and no CDM-motivated priors were considered. Rather Figure 1 and Figure 2 pertain to the specific case of the prior assumed in V max .
The above situation does, however, present a dilemma when considering priors in the quantity V max : the actual V max prior for the observable galaxy is likely to be more shallow than the prior that comes from the predicted CDM V max prior for all substructure (uniform in V −3 max ), but more shallow priors will give more statistical weight to parts of V max parameter space not well constrained by data. Because of these factors, the best that can be achieved with line-of-sight data is a lower limit to the model parameters. Thus, in this paper we use the CDM prior P(V max ) ∝ V −4 max dV max with an imposed cut-off of 3 km s −1 . This cut-off seems reasonable, as ∼ 3 km s −1 is expected to be a conservative lower bound to the V max values below which gas is able to condense into halos [80]. As we show, the imposed cut-off of 3 km s −1 does not affect the predicted flux from Draco, in which the posterior is dominated by the line-of-sight data. However, for Segue 1 (with only 24 stars), the issue is more subtle, as the CDM prior becomes more dominating with decreasing values of V max . This is would be particularly true if the V max cut-off extended down to arbitrarily low values, though our choice of a (physically-motivated) cut-off somewhat reduces the effect of the prior, even for Segue 1.
In summary, our prior for the CDM halo parameters is given by with P(r max |V max ) as in Equation (16) and P(V max ), P(1/η) according to the prescriptions given above. Finally, we have little physical basis to choose the prior for the anisotropy parameters V . Additionally, these parameters are not well constrained by the data. Fortunately, these parameters, unlike the halo parameters, do not have a direct effect on the derived flux or mass (only indirectly through the Jeans equation). Thus, we expect that the choice of these priors will not have as large of an impact on the result as compared to the halo priors. Figure 3 shows the difference between the Segue 1 mass withing 300 pc using our prior assumption and the mass assuming an isotropic velocity distribution. Although our prior assumption biases the probability distribution to lower mass values, the effect is not as extreme as with the V max prior.
The prior choices for all of the parameters (including particle physics quantities) considered in this paper are shown in Table 1.

Monte Carlo Markov Chain methodology
Here we review the MCMC formalism necessary for our analysis; we refer to the papers referenced for more details.
The goal of an MCMC algorithm is to generate a series of points in parameter space (called "a chain") with the property that their density is distributed according to the posterior pdf one wishes to explore. Then from the chain of "accepted" points, the marginal probability Outer anisotropy distribution for each of the parameters is recovered by simply binning the points in the chain, and ignoring the uninteresting coordinates (two-dimensional distributions are obtained in a similar manner). In our case, we wish to explore the joint parameter space spanned by the particle physics model parameters, C and N , and by the dwarfs model parameters, H and V . So we are dealing with a total of 15 parameters. Because the CMSSM likelihood and the stellar kinematics likelihood are independent, the joint posterior factorizes as (more details are given below) A great advantage of the MCMC procedure lies in it efficiency, for the computational effort scales roughly proportionally with the dimensionality of the parameter space being explored, rather than exponentially. The true power of MCMC methods, which we specifically utilize in this paper, lies in the fact that in addition to obtaining the distributions for the model parameters, the probability distribution for any function of the model parameters is obtained by simply determining the function at each of the accepted points in the chain. In this way, we may easily determine the distribution of any parameter that is derived from our base set of parameters by post-processing the chain of accepted points. Examples of derived parameters that we are interested in for the purposes of this paper include the dark matter mass of the dwarfs or the gamma-ray flux. While the former quantity is just a function of the parameters that describe the dark matter halo, (H , V ) the latter quantity requires us to combine the probability distributions on (C , N ) with the probability distribution from the dwarf kinematics. To understand how these probabilities can be combined, and quantities such as the flux can be robustly calculated, we appeal to the general properties of conditional distribution functions. A Markov chain of the joint posterior distribution P (C , N , H , V |D, d), Equation (18), can be obtained from the combination of the two posteriors P(C , N |D) and P(H , V |d) as long as the two joint probability distribution can be factorized as in Equation (18) above. This is, in fact, the case here since the likelihood for the particle physics model is unaffected by the stellar kinematic data. It is also the case that the halo density profile and stellar velocity dispersion anisotropy parameters are not affected by the particle physics constraints, hence the two likelihoods are independent. On the other hand, if we had included, for example, gamma-ray flux upper limits from Imaging Air Cerenkov Telescopes (ACTs) [81] in the particle physics likelihood, this would couple the two separate parameter spaces and invalidate the above decomposition. At present including such data is not necessary because the upper limits are well above the CMSSM parameter space, however it will be desirable and indeed necessary in the future.
To obtain a Markov chain on the dSph model parameters, we opted for a combination of the slice sampling [82] and the Metropolis-Hastings [83] algorithm. The advantage of slice sampling is that, unlike the Metropolis-Hastings algorithm, its efficiency is not strongly linked to the proposal pdf. Whereas, with a good proposal pdf, the Metropolis-Hastings algorithm converges faster to the desired distribution. Thus, we obtain an initial proposal pdf using slice sampling and then employed it in the Metropolis-Hastings algorithm to derive our final posterior. (For the actual slice sampling methodology, see Neal [82] and Lewis and Bridle [84]. For the Metropolis-Hastings methodology see Christensen et al. [85], Lewis and Bridle [86], Baltz et al. [87], and Lewis and Bridle [84] for details).
In practice, we found that our slice sampling run took 3-4 likelihood evaluations per point and offered fairly good convergence. However, our subsequent Metropolis-Hastings run had an acceptance rate of 30%-50% with excellent convergence. Nine chains of 30,000 points were obtained, thinned, and then combined with the SuperBayeS chains (as per the method outlined in Appendix A.2). We refer to Appendix A.1 for a more detailed discussion of the convergence criteria we use for our chains.

Boost from Halo Substructure
In the previous sections, we have outlined our modeling of the dwarf halos from stellar kinematics and our method for scanning the CMSSM parameter space and used this modeling to predict the flux under the assumption of a smooth halo. The final ingredient we must add to the flux predictions is the boost from halo substructure. The goal of this section is to derive the probability distribution for the boost factor, accounting in the most reasonable possible manner for the astrophysical and particle dark matter uncertainties that enter the calculation.

Defining Boost
Dark matter halos form hierarchically, and this results in a population of surviving gravitationally-bound substructure. High resolution dissipationless numerical simulations reveal substructures in z = 0 Milky Way-size halos with a mass spectrum that rises towards smaller masses with N(> m) ∝ m −α and α ∼ 0.9, down to the smallest masses that can be resolved, currently m 10 −6 M host ∼ 10 6 M ⊙ [27,79]. However, as we discuss below, this mass is some ten orders of magnitude larger than the minimum mass expected for CDM structure, m min . While numerical simulations that focus on nested regions within very high redshift halos have demonstrated that halos with masses close to the filtering mass do survive the initial process of halo formation [88,89], more modeling will be required to better understand the survival probability, density structure, and precise mass spectrum of the smallest CDM substructures at z = 0.
Because substructures themselves were assembled from smaller units prior to infall into their host, we expect a hierarchy of substructures within substructures that extends down to m min . The first explorations of this hierarchy of mass functions and substructure distributions is just now becoming viable in state-of-the art simulations [27]. The dark matter halos of dSphs are substructures of halos of the Milky Way, therefore their substructure fractions and mass functions are less well explored, but there are clear expectations. Relative to substructures within the Milky Way halo, substructures within dSph halos are expected to be older, and this gives them more time to assimilate into the smooth component of their hosts. Moreover, the halos of dwarf satellites do not get replenished by accreted field halos to the extent that an isolated field halo would. Both of these effects will act to reduce substructure fractions in dSph halos relative to isolated galaxy halos.
Following previous treatments, we define the boost B such that total gamma-ray flux Φ from a halo of mass M is Φ(M) = [1 + B(M, m min )]Φ(M). HereΦ is the flux that comes from a smooth halo of mass M and and the boost includes a contribution from all subhalos larger than m min , which is set by particle physics. We begin with the formulation in Strigari et al. [10], but note that the formalism of Pieri et al. [25] provides similar results. Adopting the above definition, we may write the boost as an integral that accounts for substructures going down the CDM mass hierarchy [10]: Here we have assumed that halos of mass M host substructures of mass distribution dN/d ln m = A(M/m) α for m < qM and (for now) we assume a self-similar substructure hierarchy. Written in this way, we can see that the total boost in gamma-ray luminosity depends sensitively on a competition between the smooth fluxΦ(m), which tends to decrease towards smaller masses, and dN/d ln m, which rises to small masses. For halos described by NFW density profiles [e.g. 30, 74] with a concentration-mass relation that follows c(m) ∝ m −µ with µ ∼ 0.1, we expectΦ(m) ∝ m ξ with ξ ≃ 1 − 2.2µ [10]. Setting q = 1, we find that the boost B(M, m min ) is approximately proportional to (M/m min ) α−ξ . Since the ratio (M/m min ) ∼ 10 15 for host halos of relevance, the precise value of the quantity α − ξ (∼ 2.2µ) at the smallest mass scales becomes crucial in determining whether the boost is significant or negligible. By making the most optimistic assumptions possible (q = 1, self-similar substructure hierarchies, and optimistic assumptions for the density structures of small halos with µ ∼ 0.1) Strigari et al. [10] and Kuhlen et al. [26] have shown that the boost should be no larger than about 100, and more typically 10 for m min ∼ 10 −6 M ⊙ .

Analytic Solution for Boost
The above solutions to Equation (19) [10,26] were obtained by numerical integration. Here, we present an analytic solution to Equation (19), and discuss its utility. We begin by assuming that the flux at a fixed host halo mass scales roughly as a power law,Φ(M) ∝ M ξ . With this, we may rewrite Equation (19) as  (20) can be recursively solved as: We may now solve for B n (M, m min ) by taking the Laplace transform of Equation (22) (see Appendix B). After inversion, this yields (for n > 0) where γ(a, x) is the lower incomplete gamma function defined by In the above analysis, we have assumed that the mass function is self-similar for all levels of substructure. However, this is unlikely to be the case and simulations [27,79] see less substructure in subhalos than in host halos. To account for the fact that the mass function may differ for each level of substructure, it would be useful to perform the same analysis for the case where the mass function varies independently for each level of substructure. We define the mass function to be at level i. Here, i = 0 would apply to the parent halo whereas i = 1 would apply to the subhalos and so on. Also, we let q → q i following the same notation. Using the same analysis as above, we can rewrite Equation (23) as For completeness, we present the α i = ξ solution : The utility of the analytic solution derived above stems from the fact that it lets us account for each level of substructure separately. We find that in the most extreme circumstance of ξ −α = −0.2 and m min = 10 −10 M ⊙ , inclusion of only the subhalos and sub-subhalos leads to an accuracy of 98.4%. Thus, it is not necessary to go below n = 2, for accurate boost predictions -for most part it is sufficient to just resolve subhalos to estimate the boost.

Annihilation luminosity -mass relation
In order to compute the boost, we need to estimate the gamma-ray luminosity from a subhalo of the halo under consideration. The luminosity depends strongly on the log-slope of the inner density profile of the subhalo and the effect of the tidal forces of the parent halo. In order to estimate these effects, we consider a simplified model of subhalo distribution and tidal effects.
We consider two possibilities for the spatial distribution of the subhalos. One is motivated by the current high-resolution simulations of Milky Way sized halos [75,90,91,92], where it was found that the number density of subhalos with masses greater than ∼ 10 −6 the mass of the parent halo falls off as 1/r 2 at large radii. At small radii, tidal forces of the parent halo will destroy the subhalo and hence reduce the population in the inner most regions. However, it is not clear that the smallest subhalos extending down to earth mass subhalos should have this, more diffuse, spatial distribution. We therefore also consider the other extreme where the spatial distribution of the subhalos follows the smooth distribution. We choose the second possibility for the final computations because the boost is primary affected by the smallest mass halos. We randomly generate subhalos according to the chosen spatial distributions.
Before it fell into the host halo, we allow for the log-slope of the subhalo inner density profile to range between a = 0.8 and a = 1.4, motivated [93]. We note, however, that more recent simulations indicate that there is no asymptotic inner slope [27,74,94]. The density profile assumed for the subhalo when it was outside the host halo is where a = 1 corresponds to the NFW profile and the concentration is defined by the ratio of the virial radius and the scale radius as c = R vir /r 0 . We assume a power-law virial mass function for the subhalo before it fell into the host halo and generate masses randomly for the subhalos. Since we start by picking the virial mass of the halos, we use the observed mass function for field halos. However, there is a huge extrapolation involved -we assume that the power-law extends all the way down to the minimum halo mass where the mass function gets truncated. Given the virial mass, we find the concentration of the halo using the field concentration-mass relation including the scatter. While these relations were derived for the NFW profile, we adapt them for the profile in Equation (28). Assigning a concentration is the most uncertain part of the calculation; in fact, as we will see, the flux depends sensitively on the concentration of the smallest halos.
Following the above procedure, we are able to set r 0 and ρ 0 , and from this we then find the corresponding ergodic distribution function, f a (E). Given the implied tidal radius, we adopt the following simple model for tidal effects wherein any given tidal radius defines a relative energy E 0 below which the distribution function of the subhalo drops to zero. This is analogous to the King models for the stars. That is, Two points are worth noting: (1) E 0 is fully specified by the ratio of the tidal radius to r 0 , and (2) the product ρ 0 r a 0 has to be unchanged because the tidal effects do not change the innermost regions of the subhalo. We find that the density profiles resulting from the distribution function in Equation (29) are well fit by where r t is defined here as the tidal radius and n is a function of E 0 (or equivalently r t /r 0 ). For any given ρ 0 and r 0 , this model then provides a relation between r t and the tidally truncated mass of the subhalo. If r t is smaller than r 0 , we consider the subhalo destroyed.
To fix the mass and tidal radius for a given ρ 0 , r 0 and distance from the center of the host halo, we iteratively solve the Jacobi relation between r t and mass, under the approximation that both the host halo and the subhalo are isothermal. We allow for the fact that the satellite could have been closer in its orbit by choosing the closest distance to the center along the orbit to be a factor less than one, chosen randomly from 0 to 1, times the current distance.
This is now a fully specified model and we may now predict the gamma-ray annihilation luminosity as a function of the tidally truncated mass of the subhalo. At the high mass end of the subhalo mass function, the predictions will have large scatter because of sample variance. At smaller masses, where most of the flux comes from, the large number of subhalos contributing to the flux means that the predictions essentially have no scatter.
As a means of illustrating the important uncertainties associated with subhalo density structures we explore two simple models for the concentration-mass relation. The first is a power law (PL) model c(m) ∝ m −0.1 that is normalized to match the results of N-body simulations over the mass ranges that have been probed by direct simulations (m 10 8 M ⊙ ). The second is the simple analytic model Bullock et al. [30] with the specific implementation of Maccio' et al. [95]. This model (denoted B01 below) also matches simulations down to the mass-scales explored by simulations, but links the value of c to the spectrum of CDM density fluctuations via estimated collapse times ‡. While the power-law extrapolation leads to concentrations of order 500 for earth-mass halos, the B01 model concentrations are much lower, and plateau for small masses because the slope of the power spectrum of density fluctuations is more shallow on small scales.
We note that the flux depends primarily on the profile around r 0 and smaller. The concentration sets ρ 0 since ρ 0 ∝ c 3 and the flux is proportional to ρ 2 0 r 3 0 . We find that the scaling formula suggested by Strigari et al. [10],Φ(M) ∝ M c(M) 2.2 , works well (here M is the subhalo tidally-truncated mass) for power-law c(M) functions. We did not find any systematic differences in this scaling for the two assumptions about the spatial distribution of the subhalos. However, we do find that the more concentrated spatial distribution (NFW) implies systematically larger fluxes (by a factor of about 2). In the self-similar calculation below we only include the scaling information.
We do not consider the effect of the assumptions about the spatial distribution of subhalos on the signal-to-noise. In particular, if the distribution of the subhalos is as shallow as 1/(r s + r) 2 , then much of the signal (even for more moderate boosts of order 1) will come from the "outer" regions. As we integrate the signal outwards from the stellar core, two effects must be considered. One, the signal-to-noise depends sensitively on the angular acceptance region around the satellite. If the background is mainly extragalactic (or at least constant across the angular region of the galaxy), it will scale linearly with the angular region covered. The scaling of the signal with the angular acceptance can only be estimated if the spatial distribution of the subhalos of a satellite is known. Two, we cannot accurately predict the tidal radius beyond which both the dark matter density profile of the satellite galaxy is cut-off. For the boost calculation in the MCMC exploration, we assume that the spatial distribution of the subhalos follows that of the smooth halo. We reiterate that the relevant quantity is the distribution of small subhalos that cannot yet be resolved by simulations and those subhalos are the ones that contribute dominantly to the gamma-ray flux.

Minimum Mass CDM Halo
We now have an expression for the boost as a function of (1) the host halo mass, (2) the mass function of CDM subhalos, (3) the concentration-mass relationship for CDM subhalos, (4) the minimum mass CDM halo, m min . The uncertainty in the fourth item arises from the unknown cut-off scale in the mass function of CDM substructure. As mentioned above, this cut-off scale is well below the resolution of modern day CDM numerical simulations of Milky Way like galaxies at z = 0. The smallest halo size is set by the horizon at kinetic decoupling in the early universe. The kinetic decoupling in turn depends on the scattering interactions ‡ The choice of B01 is representative of many proposed analytic models that link c to the power spectrum. Among these it has the steepest faint-end slope. of dark matter with standard model fermions, as well as the free-streaming after decoupling [33,34,35,96], Here g e f f = 10.75 is the effective degrees of freedom when the CDM particle freezes out, and where |M(t → 0)| is the matrix element of elastic scattering in the limit of no momentum transfer (t → 0). Bertschinger [35] has calculated m min assuming a Bino-like neutralino, finding typical values m min ∼ 10 −4 M ⊙ for typical WIMP masses. We follow the kinetic decoupling calculation of Bertschinger [35] but generalize it to include Wino and Higgsino-type neutralinos by determining C at each point in the SuperBayeS chain. In Appendix C, we present the relevant Feynman rules and matrix elements that enter into calculating the scattering matrix in Equation (33).
It is important to note that in our analysis of the scattering matrix element, we assume that before kinetic decoupling occurs, the WIMP interacts solely with electrons and neutrinos. In doing so, we have neglected the effects of muon scattering, which will be relevant if kinetic decoupling occurs at temperatures comparable to the muon mass. Including the nonnegligible muon abundance would in turn modify Equation (32), which is beyond the scope of our present analysis.
In Figure 4, we show the resulting posterior distribution of m min accounting for the entire presently viable parameter space of the CMSSM. We find ∼ 95% c.l. values for the minimum halo mass within the range ∼ 10 −9 −10 −6 M ⊙ . The features in the likelihood for m min arise from the probability distributions in the CMSSM parameters; due to the non-linear transformation between these parameters and the minimum mass halo in Equation (32) any small features in the CMSSM parameters are strongly exacerbated in the m min likelihood.

Boost Predictions: Two models
Running Monte Carlo simulations (see section 5.3), we found that the luminosity Φ(M) has a power law behavior for both the PL and B01 models. In terms of the first-order estimates of the boost discussed in association with Equation (19) above, the relevant combination of the mass-function log-slope (−α) and the luminosity function log-slope (ξ) takes the value α − ξ ≃ 0.2 and 0.1, for the PL and B01 models, respectively. Figure 5 shows the boost for a 10 8 M ⊙ halo across its relevant range for both models. The ensuing posterior pdf for the boost for Segue 1 for these two models is shown in Figure 6. The difference in boost between these two models is an order of magnitude and this underscores the large effect that the (uncertain) concentration of the lowest mass halos has on the overall boost. This leads to a natural uncertainty for any flux prediction and emphasizes the need to measure the concentration and mass function well for halos of the smallest size. In order to be conservative, for the remainder of this paper we assume α − ξ = 0.1.   (23)) and the B01 (dashed red) concentration models. Note that both of these concentration-mass models are consistent with current simulations.

Flux Predictions and Detection Prospects
In the analysis above we have assembled the necessary ingredients to make gamma-ray flux predictions given constraints from the halo kinematics and the CMSSM. The corresponding kinematic constraints on the (smooth) dark matter halos of Draco and Segue 1 are illustrated in Figure 7. Now, as a final step in the process of making the flux predictions, we must specify the angular region around the center of the dSph within which the flux is calculated, .i.e. θ max in Equation (1). Determining the most optimal value of θ max which maximizes the signal and minimizes the background contribution depends on the experiment and backgrounds under consideration. Satellite experiments such as Fermi [97] and ACTs, such as HESS [98] and MAGIC [99], are expected to have an angular resolution of ∼ 0.2 • and 0.1 • , respectively. If the extent of the gamma-ray emission region for the dSphs is larger than these point-spread  functions the dSphs would be resolved as point sources.
Our approach in this paper is to provide an algorithm for determining the optimal value of θ max , given the best fitting halo parameters and the measured background spectra. More specifically, we determine the optimal value of θ max from our Markov chains, along with the estimation of the detection significance, σ = N s / √ N s + N b , where N b is the number of background photons and N s is the number of signal photons within a given angular acceptance. We note that the definition of detection significance here is only meant to be an approximation; more realistically the maximum signal-to-noise will depend on the angular distribution of the signal and background, and detailed detector specifications. Nonetheless, proceeding with the above definition, at each point in our chains, we calculate the value of θ max for which the significance σ is maximized, and then construct the pdf of the quantity, θ max . In order to convert the flux to a detected number of photons, we take the total exposure, defined as the orbit-averaged effective area times the observation time, to be 3 × 10 11 cm 2 s ≃ 2000cm 2 × 5years [97].
For our input background spectrum we perform a standard analysis utilizing the EGRET diffuse background measurements at high Galactic latitude. We consider just the diffuse extragalactic background, though note that including contributions from the residual Galactic component may increase the total background by flux about an order of magnitude. The diffuse backgrounds over the relevant energy range of ∼ 1 − 100 GeV for WIMP annihilation will of course be mapped with even greater precision in the very near future with Fermi. The diffuse extragalactic background is seen to fall off according to the power law dN/dE ∼ E −2.1 [100]; we integrate this background spectrum for energies > 1 GeV, and compare to the number of photons produced by dark matter over the same energy range. The resulting pdf for θ max is shown in Figure 8. As Figure 8 shows, accounting for the diffuse backgrounds, the optimal angular extent for each galaxy is ∼ 0.2 • (the PSF of Fermi), with a gradual tail that extends to larger angles.
Given the distribution of θ max for each galaxy, in Figure 9 and Figure 10 we present our results for the flux of gamma-rays from each galaxy with energies greater than 1 GeV. In these figures, we show the gamma-ray flux versus two CMSSM parameters, the neutralino mass and the spin-independent cross section. For these figures we assume the "B01" boost model, where the minimum value for the halo mass is self-consistently calculated from the CMSSM parameter space. Because it has more line-of-sight velocities, the flux distribution for Draco is more strongly-constrained. As is seen, the predicted fluxes are similar for both galaxies, despite the fact that Segue 1 is more than a factor three times closer than Draco. The main reason that these fluxes are similar traces back to the prior used in Equation (16); given the velocity dispersion and half-light radius, Segue 1 prefers to reside in a halo with  lower V max relative to Draco. The corresponding one-dimensional flux distributions are show in Figure 11.
Using the definition of significance, σ (defined above), and an assumed exposure, we can obtain a rough estimate of detection prospects. As an example, we consider the detection  [30] and power-law halo concentration models (discussed in Section 5.5) respectively. The dashed-dotted line denotes the flux above which the signal to noise is about 3, using the estimated definition of signal-to-noise is section 6. Notice the double peaked feature in the Draco posterior. This is due to the fact that the astrophysical component of the flux is very well constrained by data causing the uncertainty in the flux to be dominated by the non-Gaussian CMSSM parameter space. The uncertainty in the flux for Segue 1, though, is dominated by the astrophysics and accordingly is very Gaussian in shape.
prospects at a signal-to-noise greater than 3 (σ = 3) in five years (same exposure as given previously). Comparing to our flux predictions in Figure 11, we see that with the standard B01 boost model, the parameter space spanned by the flux posterior of Draco lies below this estimated flux limit (vertical dot-dash line) and only a small fraction of the parameter space reaches this limit for Segue 1. On the other hand, taking the optimistic PL boost model, in which the mean boost is ∼ 20, both Segue 1 and Draco exhibit significant regions of their parameter space that are detectable given our estimated Fermi sensitivity. Specifically, for the PL boost model, we find ∼ 20% of the Draco flux parameter space is detectable with signalto-noise greater than 3 with five years of Fermi observation (and similarly ∼ 13% of the Segue 1 parameter space at the same level). We note again that the assumed priors (see below) have a large effect on the predicted flux and our general strategy in this work has been to quantify the minimum expected flux. We also note that for our PL model we assumed a spatial distribution for the subhalos that tracks the underlying host halo mass distribution. Making this spatial distribution more extended will lower the flux from the central region. Our results for the B01 concentration boost model are broadly consistent with the results of Pieri et al. [16], who also considered annihilation flux from Draco but with an energy threshold of 100 MeV.
For the PL boost model, it would be interesting to also consider the flux of gammarays from annihilations in the unresolved substructure of the Milky Way and ensure that this is not in conflict with the measured EGRET background [25]. However, such an exercise is hampered by the facts that the Galactic disk could significantly change the substructure distribution as well as abundance, and there is no systematic way of taking this into account Of course, to extract a dark matter annihilation signal above background, a detailed comparison between all of the input spectra is required, rather than simply counting photons above some energy threshold. In order to make a detailed comparison between spectra, it will be important to determine the characterstic energy for photons from dark matter annihilation. We define the characteristic energy as E max , the energy at which the quantity E 2 dN/dE peaks. For each point in our Markov chain we determine the characteristic energy, and the resulting pdf for E max is given in Figure 12.
Finally, we re-vist our above discussion of priors in the context of flux predictions. We see from Figure 11 that when the astrophysical data is well constraining the majority of the flux uncertainty comes from the particle physics parameter space. This severely reduces the impact of the astrophysical priors on the overall uncertainty of the flux. But for galaxies like Segue 1 with only 24 line-of-sight velocities, the likelihood for its parameters is much less constraining and thus astrophysical priors have a large impact on the flux posterior. In Figure 13, we show how the assumed prior affects the resulting calculation of the flux from Segue 1 and Draco. The curves in Figure 13 are for the same priors as in Figure 1. It is thus clear that, given the small number of stars from Segue 1, the assumed prior has a significant effect on the flux calculation. The CDM V max prior forces Segue 1 to reside in halos with smaller V max and this results in smaller predicted fluxes. We note that the physics of faint-end galaxy formation should change this prior significantly and would generally push the prior towards favoring larger V max values. Thus our chosen priors that have the effect of providing the minimum expected flux from Segue 1. Future data sets for this object, and all other ultra-faint satellites, will be crucial for constraining the dark matter mass and estimating the gamma-ray fluxes.
The priors on the CMSSM parameter space may also play an important role in determining the flux expectations from the dwarfs [43]. We have not undertaken a systematic study of the effect of these priors on the fluxes here. We chose a uniform prior in all the CMSSM parameters including the masses. , and ln(V max ) (blue, dotted) is assumed. As in the case of the mass (see Figure 1), flat priors in increasing negative powers of V max causes the Segue 1 posterior to be more biased than the kinematically well-constrained case of Draco. For these figures, we have set the boost factor equal to zero.

Discussion and Conclusions
Constraining the properties of the dark matter particle in indirect detection experiments will require a firm understanding of the astrophysical uncertainties that contribute to the flux. Dwarf satellites of the Milky Way are particularly interesting targets in this regards because they are the most dark matter dominated objects known and they are largely free from astrophysical uncertainties which result from the presence of baryonic physics. In this paper, we have taken an step towards quantifying uncertainties in flux predictions from dSphs by providing a framework within which both particle physics and astrophysics uncertainties can be included at once.
We have combined our MCMC method for determining the dark matter distributions of dwarf satellites using stellar kinematics with the SuperBayeS MCMC package which determines the preferred ranges for parameters of the CMSSM. We have focused on two specific dSphs, Segue 1 and Draco, as example cases. Our methods allow us to provide a broad outline for the prospects of detection of these satellites with gamma-ray experiments such as Fermi.
The main results of our paper can be summarized as follows: • We find that both Draco (at 80 kpc) and Segue 1 (at 23 kpc) are expected to have grossly similar fluxes, though the flux from Segue 1 is subject to larger uncertainties because of its relative lack of kinematic data and smaller stellar extent. For the most conservative assumptions, Segue 1 prefers to reside in a halo with lower V max relative to Draco, and therefore it has a similar overall flux despite its relative proximity. However, for a flat prior in log(V max ), the flux from Segue 1 can be much larger than Draco (see Figure 13). This result motivates future observations of stellar kinematic data for Segue 1. We note that our results for the flux, unless otherwise stated, are based on the more conservative prior and conservative boost model, and hence they quantify the minimum expected flux for the assumed CMSSM priors.
• We have provided the first self-consistent calculation for the boost in the flux signal from halo substructure, taking into account both the CMSSM model and the recent results from high resolution numerical simulations. We show that the dominant uncertainty in the boost calculation comes from the assumed halo concentration versus mass relation for halo substructure on mass scales down to the scale of the minimum mass halo. If we assume a model that links halo concentrations to the power spectrum, then we obtain typical boost factors of order unity. If we instead assume a power law continuation of the concentration-mass relation down to the minimum mass halo, the average boost factor increases by an order of magnitude to ∼ 20 (see Figure 6). This boost would be reduced if the spatial distribution of the smallest subhalos is more diffuse that the smooth dark matter halo component, which is what we have assumed in this paper. We also note that our analytic solution for the boost shows that resolving subhalos, and in some cases subsubhalos, is sufficient to get an accurate estimate of the boost. These facts motivate future high resolution simulations of halo substructure that can measure the concentrations of small dark matter halos directly, as well as map their spatial distribution within the host.
• We have provided a broad outline of the prospects for detection of these satellites with gamma-ray experiments, focusing specifically on Fermi. We find that, given the diffuse backgrounds, the most optimal solid angle to view these galaxies is ∼ 0.2 − 0.3 • (see Figure 8). Optimistic fluxes for these galaxies are approximately a few times 10 −11 cm −2 s −1 . For the boost model with the power-law extrapolation of the concentrationmass relation and a subhalo spatial distribution that tracks the underlying host mass distribution, we estimate a ∼ 20% chance for a ∼ 3σ dark matter gamma ray signal from Draco after 5 years of observation. This expectation is highlighted in Figure 11. We note that given the observed uniformity in the central density of the dark matter halos of the dwarfs [22], it should be possible to stack them and increase the signal-to-noise.
• We have provided an updated calculation for the minimum halo mass in CMSSM. Our results broadly agree with those of Profumo et al. [34], and are typically lower than the masses in Bertschinger [35] who only considered bino-like neutralinos. We find that the minimum mass CDM halo lies in the range 10 −9 − 10 −6 M ⊙ (see Figure 6). The inclusion of this effect results in slightly larger boost estimates.
The methods presented here provide a concrete methodology for addressing the uncertainties inherent in dark matter indirect detection from both particle physics and astrophysics perspective. The focus on Segue 1 and Draco was meant to illustrate how two vastly different kinematic data affect the predictions for the gamma-ray flux from dark matter annihilations. It is very interesting to note that Fermi observations of the dwarf satellites of the Milky Way could be relevant for constraining the Supersymmetric parameter space.
where the power spectrum is normalized to a variance of one and Equation (A.1) is fitted over the range 0 < k < 10k ⋆ (see Dunkley et al. [101] for details). It should be noted however that this test does not guarantee that the chain has sampled the entire parameter space, but only that its autocorrelation lengths is sufficiently short to ensure independence among the samples.

Appendix A.2. Combination Methodology
If two parameter sets are completely independent, then the statistics of the combined parameter set is independent of whether each parameter set is derived independently or jointly. Thus, if a posterior is separable (P(φ, ψ) = P φ (φ)P ψ (ψ)), then the total posterior may be found by the combination of the individual distributions. But care must be taken combining chains with different weights. If the resulting combined weights differ from the individual distribution weights, then the resulting combined distribution runs the risk of differing from the original, pre-combined, distributions when marginalized over. A solution to this problem is to explicitly repeat points in each chain (thereby separating the weights out among the repeated points) in such a way that the two chains sets of weights coinside.

Appendix B. Boost Derivation
In this appendix, we derive Equation (23). We start by taking the Laplace transform of Equation (22): where F n (s) is the Laplace transform of D n (b). By recursively substituting F n (s), we can derive

Appendix C. Scattering Matrix Calculation
In this appendix we derive the cross section for neutralino-lepton elastic scattering. We closely follow Chen et al. [102]'s procedure in finding the neutralino-neutrino scattering matrix element with the exception that we generalized their result for any lepton. We start with the amplitudes for s, u, and t channels represented in Figure C1 (see Jungman et al. [36], Haber and Kane [103], Gates and Kowalski [104] for the appropriate Feynman rules).
(C.7) and the coupling constants X ′ ℓi jn , W ′ ℓi jn , and O ′′ nn{R,L} are defined in reference [36]. Here, the indices ℓ, i, j, and n represent the family type (charged lepton or neutrino), the lepton flavor, the slepton flavor, and the scattering neutralino. Here, C = γ 2 γ 0 is the charge conjugation matrix. Using the above amplitudes, the combined matrix element is Here we used the simplifying notation (C. 17) In the limit of highly relativistic leptons and by assuming m χ ≫ E ℓ , the specific momentum scalar products can be approximated by This leads to a matrix element of