Abstract

Both major galaxies in the Local Group host planar distributions of co-orbiting satellite galaxies, the Vast Polar Structure (VPOS) of the Milky Way and the Great Plane of Andromeda (GPoA). The Λ cold dark matter (ΛCDM) cosmological model did not predict these features. However, according to three recent studies the properties of the GPoA and the flattening of the VPOS are common features among subhalo based ΛCDM satellite systems, and the GPoA can be naturally explained by satellites being acquired along cold gas streams. We point out some methodological issues in these studies: either the selection of model satellites is different from that of the observed ones, or an incomplete set of observational constraints has been considered, or the observed satellite distribution is inconsistent with basic assumptions. Once these issues have been addressed, the conclusions are different: features like the VPOS and GPoA are very rare (each with probability ≲10−3, and combined probability <10−5) if satellites are selected from a ΛCDM simulation combined with semi-analytic modelling, and accretion along cold streams is no natural explanation of the GPoA. The origin of planar dwarf galaxy structures remains unexplained in the standard paradigm of galaxy formation.

1 INTRODUCTION

In recent years, the notion that satellite objects of the Milky Way (MW) are preferentially distributed and orbit within a common plane (Kunkel & Demers 1976; Lynden-Bell 1976), a vast polar structure (VPOS),1 has been substantiated by additional and improved observational information. The distribution of the MW satellite galaxies can be described by a relatively thin (about 20 kpc root-mean-square, rms, height) plane (Kroupa et al. 2005; Metz, Kroupa & Jerjen 2007; Pawlowski, Kroupa & Jerjen 2013), whilst their proper motions indicate that up to 8 out of the 11 brightest (termed the ‘classical’ ones by Metz et al. 2008) MW satellites orbit in the same sense (‘co-orbit’) within the structure (Metz et al. 2008; Pawlowski & Kroupa 2013). Globular clusters and stellar and gaseous streams appear to preferentially align with the VPOS too (Forbes et al. 2009; Keller, Mackey & Da Costa 2012; Pawlowski, Pflamm-Altenburg & Kroupa 2012a).

A similar planar distribution of satellites has been sought around M31 (Koch & Grebel 2006; McConnachie & Irwin 2006; Metz et al. 2007), but it required the homogeneous data set of M31 satellite distances provided by Conn et al. (2012) to allow Ibata et al. (2013) and Conn et al. (2013) to discover such a thin (of about 10 kpc rms height) planar satellite structure with high statistical significance. Like the VPOS, this Great Plane of Andromeda (GPoA)2 appears to be co-orbiting, but it consists of only about half the currently known M31 satellites. Similar alignments of dwarf galaxies have since been identified among the secluded dwarf galaxies in the Local Group (LG) (Bellazzini et al. 2013; Pawlowski et al. 2013; Pawlowski & McGaugh 2014) and aligned systems of satellites and stellar streams are also being discovered around more distant galaxies (Galianni et al. 2010; Duc et al. 2011; Paudel et al. 2013; Karachentsev et al. 2014).

The existence of the thin, planar VPOS is a challenge for the standard Λ cold dark matter (ΛCDM) cosmological model. Kroupa et al. (2005) noted for the first time that the distribution of dark matter subhaloes in cosmological simulations is isotropic to first order, thus being incompatible with the planar VPOS. This realization found further support in the fact that most of the luminous MW satellites orbit in the same direction within the VPOS (Metz et al. 2008; Pawlowski & Kroupa 2013).

The phase-space correlation of both satellite galaxy structures (VPOS and GPoA) has led to the suggestion (Metz & Kroupa 2007) that a large fraction of the observed satellites of the MW and of M31 are in fact tidal dwarf galaxies (TDGs; e.g. Barnes & Hernquist 1992; Bournaud, Duc & Emsellem 2008) formed in the debris shed by a past major galaxy encounter rather than being dark matter dominated primordial dwarf galaxies. The formation of TDGs out of the debris of galaxy collisions has been observed to take place even in the current Universe and it can be considered observationally and numerically confirmed that this process results in phase-space correlated satellite systems (e.g. Bournaud et al. 2007, 2008; Wetzstein, Naab & Burkert 2007; Duc et al. 2011, 2014). TDGs have been found to be indistinguishable from dwarf elliptical galaxies (Dabringhausen & Kroupa 2013), are able to survive their initial burst of star formation (Recchi et al. 2007; Ploeckinger et al. 2014), and might also be responsible for the larger scale dwarf galaxy planes in the LG (Bellazzini et al. 2013; Pawlowski et al. 2013; Pawlowski & McGaugh 2014).

A number of possible galaxy encounter scenarios in the LG have been suggested to explain the formation of the observed dwarf galaxy structures, such as the tidal disruption of a Large Magellanic Cloud progenitor (Kunkel & Demers 1976; Lynden-Bell 1976; Pawlowski, Kroupa & de Boer 2011) and a past encounter between the MW and M31 (Sawa & Fujimoto 2005; Pawlowski et al. 2012a; Zhao et al. 2013). The currently most advanced scenario of such an event is that of a major merger in M31 by Hammer et al. (2010). Their numerical models were initially tailored to reproduce features observed in M31, but also resulted in the formation of TDGs from the tidal debris. Some of these turned out to reproduce the later observed GPoA (Hammer et al. 2013), while others are flung out towards the MW where they form a VPOS-like structure (Yang & Hammer 2010; Fouquet et al. 2012), which is consistent with the alignment of the GPoA plane towards the MW.

The virtue of the suggested TDG origin is that, for a given scenario, it is possible to make detailed predictions of the orbits of TDGs formed in such an encounter. This makes the different scenarios of TDG formation testable, in contrast to the scenario in which the satellites have been accreted as primordial dwarf galaxies in an essentially stochastic manner, which results in an inherent inability to predict their exact orbits. However, since being identified as a potential problem for the ΛCDM cosmological model, numerous attempts have been made to reconcile the observed satellite structures with a cosmological origin.

A proposed scenario is the accretion of groups of satellite galaxies, which would naturally result in them having similar orbital parameters (D'Onghia & Lake 2008; Li & Helmi 2008). However, Metz et al. (2009) have shown that this explanation is inconsistent with observed properties of dwarf galaxy groups. They are too extended to form narrow structures such as the VPOS or the GPoA, both of which have a thickness of ≲50 kpc.

Lovell et al. (2011) claimed that quasi-planar distributions of coherently rotating satellites similar to the VPOS arise naturally in simulations of a ΛCDM universe. However, using their data Pawlowski et al. (2012b) tested how often the orbital poles (directions of angular momenta) of subhaloes are as concentrated and as perpendicular to the Galactic disc normal as those of the observed MW satellites. The authors were unable to find support for the conclusion by Lovell et al. (2011). Instead, an orbital pole distribution similar to the observed one is extremely unlikely if the satellite orbital poles are drawn from those of cosmological subhaloes, but very likely if the orbital poles are drawn from those of tidal debris, as would be the case for TDG satellites.

Libeskind et al. (2009) and Deason et al. (2011) are two of the few studies that in addition to using the spatial information of dark matter subhaloes from cosmological simulations also compare the orbital poles of subhaloes with those of the observed satellite galaxies. Libeskind et al. (2009) find that in approximately 35 per cent of the simulated satellite galaxy systems they consider, the orbital poles of at least 3 out of 11 satellite galaxies point along the short axis of the galaxy distribution. Deason et al. (2011), searching for alignments of 3 out of only 10 satellite galaxies, find this to be the case for approximately 20 per cent of their simulated systems. However, as we have shown in Pawlowski & Kroupa (2013), additional and improved proper motion data for the 11 brightest MW satellites result in a much more concentrated distribution of orbital poles for the co-orbiting satellites, such that now at least six of the orbital poles are well aligned. This is an extremely rare finding for ΛCDM-based satellites according to fig. 9 of Libeskind et al. (2009) and fig. 8 of Deason et al. (2011).

Three new attempts at reconciling the observed GPoA and VPOS with a purely cosmological origin have been made in the past two years. These studies claim that: (1) structures such as the GPoA in M31 are not rare among satellite systems in the Millennium II simulation and thus the GPoA is consistent with ΛCDM (Bahl & Baumgardt 2014, BB14), (2) structures such as the VPOS in the MW are in agreement with a ΛCDM universe because flattened satellite distributions as the one around the MW are frequent in the Millennium II simulation (Wang et al. 2013, W13), (3) the GPoA can be naturally explained if host galaxies acquire not only their gas but also their satellites along cold streams and that the GPoA can therefore be interpreted as indirect observational evidence of the paradigm of cold mode accretion streams (Goerdt & Burkert 2013, GB13).

When searching for GPoA or VPOS analogues in cosmological simulations, it is important to decide which are their fundamental parameters. The most obvious ones describe the overall spatial distribution. For example the rms height of the satellites around their best-fitting plane (rper) and their radial distribution in that plane measured by the rms of their distances projected into the plane (rpar). An alternative is using the relative shape measured by the rms axis ratios (c/a and b/a). But the satellite planes are structures in six-dimensional phase space. The VPOS is rotationally stabilized (Metz et al. 2007; Fouquet et al. 2012; Pawlowski & Kroupa 2013) and the coherent velocity directions of the satellites belonging to the GPoA were, in fact, the most important part of the discovery of the structure (Ibata et al. 2013). When comparing simulations to observations, an analysis would therefore be incomplete if it would focus on only part of the observable dimensions such as part of the position space by only considering rper, in particular when velocities have been used observationally to detect the structures (e.g. Metz et al. 2008; Libeskind et al. 2009; Deason et al. 2011; Pawlowski et al. 2012b; Ibata et al. 2013; Pawlowski & Kroupa 2013) and are available from the simulation.

The minimum required parameters can thus be summarized by trying to give the shortest possible characterization of our current knowledge of the two satellite structures in the LG: (1) the GPoA is thin (rper) but extended (rpar) and has coherent velocity directions (ncoorb, number of satellites showing the same velocity trend). (2) The VPOS is thin (rper) but extended (rpar) and the satellite orbital poles are strongly concentrated (Δstd, spherical standard distance of the concentrated orbital poles) in the direction of the plane normal (θVPOS, angle between plane normal direction and average direction of concentrated orbital poles). In order to be as conservative as possible these parameters will only be considered as limits in the following. For example, the rms height of a simulated satellite plane is not required to match the observed values closely, but only to be smaller or equal to the observed value.

This paper is structured as follows. In Section 2 we search for GPoA-like planes in the same simulation data set used by BB14. We use the same simulation data set in Section 3 to search for VPOS-like satellite planes, as in the study by W13. In Section 4 we test whether cold accretion streams can naturally explain the GPoA, as claimed by GB13. We close the paper with concluding remarks in Section 5.

Sections 2 to 4 are all similarly structured. They start with a discussion of the previously applied method and a comparison with the available observational data, followed by the description of the method applied by us. The results of this new method are presented, ending in a summarizing discussion. For observations, if not mentioned otherwise, in the following we will use the positions of M31 and its satellite galaxies based on Conn et al. (2012), while all remaining positions and the line-of-sight velocities are based on the catalogue of LG galaxies by McConnachie (2012). All properties of the simulated satellite systems correspond to the situation at a redshift of z = 0.

2 SEARCHING FOR A M31 GPoA-LIKE STRUCTURE IN ΛCDM SIMULATIONS

2.1 Previous method and its limitations

As mentioned in Conn et al. (2013), a satellite galaxy plane such as the GPoA is not expected in a ΛCDM cosmology. The discussion in Ibata et al. (2013) also illustrates that there is no obvious explanation for this structure. However, both studies have focused on the observed satellite distribution only. How likely a structure such as the GPoA is in the standard galaxy formation model based on the ΛCDM cosmology can be quantified by comparing the observed satellite system with simulations of M31-like galaxies. Herein lies the importance of the approach taken by BB14. They have searched for GPoA-like structures among simulated satellite galaxy systems based on the Millennium II simulation (Boylan-Kolchin et al. 2009). This dark-matter-only simulation has been scaled to match the cosmological parameters of Wilkinson Microwave Anisotropy Probe 7 (WMAP7; Jarosik et al. 2011) and has then been populated with galaxies using the semi-analytic galaxy formation models by Guo et al. (2013). For details we refer the reader to these papers.

BB14 conclude that the GPoA is not in conflict with the standard cosmological framework. According to their analysis, satellite configurations as thin as the observed GPoA are quite common around galaxies with a mass comparable to M31 and a significant fraction of the galaxies in the simulated satellite systems are co-orbiting.

A statistical test has to compare consistent samples. This means that the selection applied to the observed satellite system also has to be applied to the simulated one. This might only be possible approximately, but the selection should be made as similar as possible. In addition, it has to be made sure that the selection applied to the simulated satellite system is also applied to the observed satellites before the plane parameters are determined for both, preferentially by using the same method. Only then can a comparison of the plane parameters yield meaningful results. The analysis of BB14 is, in this respect, limited. We have identified the following problems.

  • Three criteria were introduced which must be met by a satellite system to be similar to the observed GPoA, based on the properties of the best-fitting (GPoA-like) plane consisting of a subset of the satellite galaxies (15 out of 27 satellites in BB14). These properties are (1) the rms height rper of the 15 closest satellites from the best-fitting plane, (2) the rms distance from the centre of the host galaxy, rpar, of these 15 satellite galaxies projected into the best-fitting plane and (3) the number of co-orbiting satellites (ncoorb). A plane in the simulated satellite systems has to be at least as narrow as the observed one. The plane also has to be at least as radially extended as observed (rparsimrparobs), otherwise rper might be small because the whole satellite system is very radially concentrated but not necessarily flat. Finally there have to be at least as many co-orbiting satellite galaxies as observed in the GPoA. Despite introducing these three independent criteria, BB14 discuss combinations of only two of them (either rper and ncoorb or rper and rpar). As only a small fraction (2 per cent) of the satellite systems fulfils the two-criteria combinations3 in the first place, it is important to determine how much lower this fraction drops for systems reproducing all three properties of the observed GPoA.

  • The Pan-Andromeda Archeological Survey (PAndAS) area was approximated as a cylinder with a radius of 128 kpc (Bahl, private communication). Drawing the simulated satellites from a cylindrical volume assumes an infinite distance to the MW, in contrast to the real PAndAS survey which is better approximated by a conical volume with the MW in the apex.

  • The used method does not guarantee that they only identify (oblate) planes, because narrow prolate distributions of satellites can have small rper too. In addition, they apparently consider a too-small PAndAS footprint of 128 kpc radius at the distance of M31, while Ibata et al. (2013) and Conn et al. (2013) clearly state that the radius is approximately 150 kpc. This results in a 27 per cent smaller area from which simulated satellites are selected, which makes the overall satellite distribution more narrow. As prolate distributions are also accepted as planes, this method is biased towards finding narrower structures (having smaller rper) which artificially increases the agreement with the observed situation.

  • 27 satellites are drawn from an area around M31 that is supposed to model the PAndAs survey footprint. However, while the major part of the PAndAS survey indeed covers an approximately circular area of 150 kpc radius at the distance of M31, the survey also covers an area of approximately 50 kpc radius around M33. This is important as two of the 27 M31 satellites covered by PAndAs are outside the 150 kpc zone: M33 and Andromeda XXII. Both objects are not associated with the GPoA. An alternative to using a circular footprint would have been to use the exact PAndAs footprint from different random viewing angles (Ibata et al. 2014). However, if one chooses to stay as close as possible to the original method of BB14, i.e. simulated satellites are only drawn from within the 150 kpc zone, then only 25 satellites must be drawn from the simulation. Therefore, instead of searching for planar alignments of 15 out of 27 satellites, BB14 would have had to search for alignments of 15 out of 25. Using a too-large number of simulated satellites increases the probability to find narrow alignments of satellites because more objects are available to select from in the same volume. As an example, the number of possible combinations of 15 satellites out of 27 is |$\frac{27!}{15! 12!} \approx 17.4 \times 10^{6}$|⁠, but only |$\frac{25!}{15! 10!} \approx 3.3 \times 10^{6}$| for 15 out of 25 satellites. Additional satellites therefore substantially increase the chance of finding a more narrow plane if the number of satellites assigned to the plane is kept fixed. The choice of the satellite number by BB14 therefore biases the results towards finding a better agreement with the observed situation.

  • Simulated satellite galaxies were selected from cubical volumes extracted from the Millennium II simulation data which are centred on the host galaxies (Bahl, private communication). These cubes have an edge length of 600/h kpc, where h = 0.704. When searching for anisotropies in the simulated satellite systems, it is crucial to draw them from a volume that is non-directional, i.e. a sphere. Otherwise, the 3D geometry of the volume from which satellites are drawn will introduce anisotropies in the data. The maximum spherical volume fitting into these cubes has a radius of rsphere = 426 kpc. However, according to figs 3 and 7 in BB14, some of their simulated satellites have distances of 700 kpc or more from their hosts. This indicates that the satellites have not been drawn from a spherical volume, but from the cubes directly, because in that case the random viewing direction will in some cases be along the corners of the cubes, resulting in a maximum distance from the host galaxy of |$d_{\mathrm{max}} = \scriptstyle\sqrt{3 r_{\mathrm{sphere}}^2}$| = 738 kpc. As the more distant galaxies are only present in the directions of the eight corners of the considered simulation cube, the inclusion of these galaxies artificially introduces a preferred directional alignment to the satellite galaxy distribution. It also results in different radial distributions, because the maximum distance of a satellite from the host depends on the orientation in which the cube is seen. This has to be avoided when trying to determine the chance of finding satellite alignments in a meaningful way.

  • As pointed out before, for the simulation data used by BB14 only satellite galaxies within a maximum distance of dmax = 426 kpc can be considered in the analysis, because no spatially unbiased data are present for the more distant simulated galaxies. However, Andromeda XVIII is at a distance of ≈460 kpc from M31 and therefore outside of the volume considered around each of the simulated host galaxies. Consequently, this galaxy must be excluded from the sample of M31 galaxies and the number of satellites drawn from the simulations also has to be reduced by one. As Andromeda XVIII is not part of the GPoA, the analysis still has to look for planes made up of 15 satellites. However, the M31 satellite galaxy Andromeda XXVII introduces a further complication here. As seen from the MW, it is superimposed on to a stream of stars. This makes a distance determination to the object difficult, such that the discovery paper only reports a lower limit of its distance, placing it at ≈760 kpc from the MW or ≈90 kpc from M31 (Richardson et al. 2011). In contrast to that, the analysis by Conn et al. (2012) places its most-likely heliocentric distance at an upper limit of ≈1250 or ≈480 kpc from M31, but its distance posteriori distribution (see their fig. 12) is strongly bimodal, indicating that the object is either at ≈1250 kpc heliocentric distance or at ≈800 kpc. If Andromeda XXVII is near, it falls within the volume extracted from the simulation, if it is far it also has to be excluded from the sample. In summary, to be compatible to the observed distribution of M31 satellites as used in the analysis of Ibata et al. (2013), one is only allowed to draw 23 satellite galaxies from the simulation data used by BB14 and to search for planar alignments of 14 of these satellites. Assuming a closer distance of Andromeda XXVII increases the numbers by one: 15 satellites out of 24. In both cases the properties of the GPoA plane fit need to be re-determined as either the satellite sample or the position of one of the satellites is different from that used in the analysis of Ibata et al. (2013).

    In summary, the two major issues are that BB14 did not demand all three of the essential GPoA criteria as they were defined in their analysis to be fulfilled simultaneously, and that their selection of model satellites differs from the selection of observed satellites in number, selection area and distance from M31. The importance of the second point can be illustrated as follows: even if there would be a simulated satellite system in which the 27 satellites with the largest baryonic mass have the exact same positions and velocities as the satellites observed around M31, the analysis of BB14 would not have been able to select all these satellites. The analysis would thus have determined different parameters from a different satellite sample than done observationally by Ibata et al. (2013), even though the distribution of the top 27 satellites would have been identical.

    Finally, we remind the reader of the following other, possibly important aspects which have not been mentioned in BB14, and which we do not consider either in our new method below.

  • Their host galaxy sample includes galaxies in high density environments like galaxy clusters, very different environments than where M31 is found. In the future, it might therefore be worthwhile to investigate whether the environment has any influence on the occurrence of co-orbiting satellite galaxy planes.

  • The M31 satellite distribution is lopsided, most of the satellite in the GPoA are on the MW side of M31 (Conn et al. 2013).

  • The GPoA is orientated almost edge-on towards the MW, being inclined by ≤5° towards our line of sight (Conn et al. 2013; Ibata et al. 2013; Pawlowski et al. 2013).

2.2 New method

For best comparability, we use the same data extracted from the Millennium II data base (Lemson & Virgo Consortium 2006) as was used by BB14. It consists of all galaxies within cubes of edge length dedge = 600/h kpc = 852 kpc centred on the main halo in the virial mass range of 1.1 to 1.7 × 1012 M. There are 1825 such main haloes. Following the recipe in BB14 we remove all main haloes that contain a galaxy with a baryonic mass larger than 7 × 1010 M within 500 kpc. Note that the simulation boxes do not reach out to 500 kpc, such that we can only be sure to remove galaxies within a sphere of radius dedge/2. We then try to reproduce their results, whilst also addressing the caveats mentioned above.4 Table 1 compares our assumptions with those of BB14.

Table 1.

Comparison of assumptions.

Parameter BB14 This work
PAndAS shape Cylinder Cone
PAndAS radius at host 128 kpc 150 kpc
Simulation volume shape Cube Sphere
Simulation volume size 852 kpc edge length 852 kpc diameter
Andromeda XXVII Far Near or far
Nsat 27 24 or 23
NGPoA 15 15 or 14
Parameter BB14 This work
PAndAS shape Cylinder Cone
PAndAS radius at host 128 kpc 150 kpc
Simulation volume shape Cube Sphere
Simulation volume size 852 kpc edge length 852 kpc diameter
Andromeda XXVII Far Near or far
Nsat 27 24 or 23
NGPoA 15 15 or 14

Assumptions of the analysis presented by BB14 (second column) and the ones in this work (third column).

Table 1.

Comparison of assumptions.

Parameter BB14 This work
PAndAS shape Cylinder Cone
PAndAS radius at host 128 kpc 150 kpc
Simulation volume shape Cube Sphere
Simulation volume size 852 kpc edge length 852 kpc diameter
Andromeda XXVII Far Near or far
Nsat 27 24 or 23
NGPoA 15 15 or 14
Parameter BB14 This work
PAndAS shape Cylinder Cone
PAndAS radius at host 128 kpc 150 kpc
Simulation volume shape Cube Sphere
Simulation volume size 852 kpc edge length 852 kpc diameter
Andromeda XXVII Far Near or far
Nsat 27 24 or 23
NGPoA 15 15 or 14

Assumptions of the analysis presented by BB14 (second column) and the ones in this work (third column).

For each host galaxy we construct a sample of satellite galaxies. We only consider satellites with masses above 2.8 × 104 M, again following the cuts by BB14. In the semi-analytic galaxy formation model which was applied to the Millennium II simulation, galaxies within subhaloes that drop below the resolution limit are allowed to survive. The position of these ‘orphan’ galaxies is determined by following that particle from the original subhalo which was most bound at the last time the subhalo was resolved. The positions and velocities of these objects are therefore less reliable than those of galaxies embedded in resolved subhaloes (Pawlowski et al. 2013; W13), but for simplicity and comparability to BB14, as well as W13 (see next section) we do not discriminate between orphan and non-orphan galaxies.

For each halo a random viewing direction is determined. Then the satellite galaxies with the largest baryonic masses (i.e. the brightest ones) are chosen following four different sample prescriptions.

  • ‘BB14-like’ sample. Nsat = 27 satellites within a cylinder of 128 kpc radius and outside of a cylinder of 32 kpc radius centred on the host galaxy are drawn from the whole simulation cube. There is no radial cut applied and thus satellites can have radii larger than dedge/2 if the viewing direction is along the corners of the cube. This sample reproduces the sample by BB14, and is included for comparison only.

  • ‘Near’ sample. Nsat = 24 satellite galaxies are drawn from a cone with the apex at a distance of 780 kpc (the distance at which we observe M31; Conn et al. 2012). The cone has an opening angle such that a circle of 150 kpc radius is covered at the position of the host (approximately resembling the main area of the PAndAS survey). No satellites closer than 2|$_{.}^{\circ}$|5 from M31 as seen from the artificial MW position are considered, in line with the procedure by Ibata et al. (2013). We apply a radial cut such that only galaxies within a sphere of dedge/2 = 426 kpc around the host can be selected. This sample is going to be compared to the GPoA parameters where Andromeda XXVII is assumed to be close to M31. We remind the reader that only 24 instead of 27 satellites are considered because the PAndAS survey area around M33, which contains two satellites in the observed sample, is not covered by the simulations and because Andromeda XVIII is too distant from M31 to fall into the spherical volume.

  • ‘Far’ sample. The third sample is similar to the second, but here only Nsat = 23 satellite galaxies are drawn. This assumes that Andromeda XXVII is far from M31, such that it is outside of the spherical volume covered by the simulated satellites. Consequently, when comparing with the observational data only planes consisting of NGPoA = 14 satellites are searched in this sample.

  • ‘Random’ sample. The last sample is identical to the third, but instead of drawing satellite positions and velocities from the simulation we first randomize both their position and their velocity direction by rotating them into a random direction. This results in them having an isotropic distribution around their hosts, but the same distribution of radial distances. Comparing the results for this and the previous two samples therefore allows us to determine the influence of the anisotropy present in the ΛCDM simulation on the results.

Note that all four samples use the same viewing direction for each halo, enabling maximal comparability. If the number of satellite galaxies found for one of the samples is too small, a new viewing direction is chosen for the halo. This step is repeated up to 10 times. If by then no viewing direction produced a sufficient number of satellites to lie within the artificial survey areas the halo is discarded from the analysis. For each halo 10 different samples for each sample prescription are chosen, utilizing different random viewing directions. Applying all beforementioned cuts, we arrive at a total of 15 000 samples of satellite systems.

Each sample is then searched for the thinnest plane consisting of NGPoA = 15 or 14 satellite galaxies. We follow a method similar to the one used by Conn et al. (2013) on the M31 satellite galaxy distribution. We determine the rms height of the NGPoA closest satellite galaxies to 10 000 planes with normal vectors that are approximately evenly distributed on the sphere, while Conn et al. (2013) and BB14 (private communication) have used random plane orientations. The galaxy sample producing the smallest rms height is identified as the thinnest plane, and a full plane fit to its NGPoA members is performed using the method described in Pawlowski et al. (2013). In short, the principal axes of the satellite system are determined by calculating the eigenvectors of the tensor of inertia for the distribution of satellite galaxy positions. The rms of the height rper of the NGPoA satellites from the best-fitting plane and their rms distance rpar after being projected into the best-fitting plane are determined. In addition, we measure the axis ratios c/a and b/a of short-to-long and intermediate-to-long axis. These quantities allow us to discriminate between elongated prolate (both ratios are small) and plane-like oblate (c/a is small but b/a large) satellite distributions, in contrast to only using the parameters rper and rpar.

In addition, we determine the number of co-orbiting satellites within the plane sample by projecting the angular momenta of the satellite galaxies on to one normal vector of the best-fitting plane. The co-orbiting number Ncoorb is then either the number of positive or negative angular momenta, whichever is larger. Note that according to this definition two satellites on very different (e.g. perpendicular) orbits are still considered co-orbiting if they both have angular momentum directions closer than 90° from the plane normal vector. The method therefore does not require the GPoA to be rotationally stabilized (which is the case if the angular momentum axes are aligned such that the structure is long lived). If measurements of the proper motions of the M31 satellites reveal that the GPoA is indeed rotationally stabilized, such that the angular momentum directions cluster close to the plane normal direction (as in the case of the VPOS), then the significance of the satellite structure would be greatly increased.

The parameters determined for each sample of simulated satellites are then compared to those determined for that observed M31 satellite sample which fulfils the same selection criteria. The parameters of the observed satellite distribution have been determined using the same routine applied to the simulated samples, in particular the same cut on the distance from M31 was applied to the observed satellite sample. We use the most likely distance moduli as homogeneously measured for all M31 satellites in the PAndAS field by means of the tip of the red-giant-branch method (Conn et al. 2012). The resulting parameters are compiled in Table 2.

Table 2.

Parameters of the observed GPoA.

Symbol Description ‘Near’ sample ‘Far’ sample
|$r_\mathrm{M31}^\mathrm{max}$| Maximum allowed satellite radius from M31 300/h kpc 300/h kpc
dAndXXVII Assumed heliocentric distance to Andromeda XXVII ≈800 kpc ≈1250 kpc
Nsat Number of satellites in sample 24 23
NGPoA Number of satellites making up GPoA 15 14
rper rms height of GPoA satellites from best-fitting plane 10.9 kpc 10.9 kpc
rpar rms distance of GPoA satellites projected into best-fitting plane 147.2 kpc 151.0 kpc
rparall rms distance of all Nsat satellites projected into best-fitting plane 136.8 kpc 138.8 kpc
Ncoorb Number of satellites in GPoA which co-orbit 13 13
c/a Minor-to-major axis rms height ratio 0.105 0.102
b/a Intermediate-to-major axis rms height ratio 0.788 0.740
Nlopside Number of GPoA satellites on MW side of M31 13 13
αMW Angle between GPoA and line of sight 4|$_{.}^{\circ}$|3 4|$_{.}^{\circ}$|9
Symbol Description ‘Near’ sample ‘Far’ sample
|$r_\mathrm{M31}^\mathrm{max}$| Maximum allowed satellite radius from M31 300/h kpc 300/h kpc
dAndXXVII Assumed heliocentric distance to Andromeda XXVII ≈800 kpc ≈1250 kpc
Nsat Number of satellites in sample 24 23
NGPoA Number of satellites making up GPoA 15 14
rper rms height of GPoA satellites from best-fitting plane 10.9 kpc 10.9 kpc
rpar rms distance of GPoA satellites projected into best-fitting plane 147.2 kpc 151.0 kpc
rparall rms distance of all Nsat satellites projected into best-fitting plane 136.8 kpc 138.8 kpc
Ncoorb Number of satellites in GPoA which co-orbit 13 13
c/a Minor-to-major axis rms height ratio 0.105 0.102
b/a Intermediate-to-major axis rms height ratio 0.788 0.740
Nlopside Number of GPoA satellites on MW side of M31 13 13
αMW Angle between GPoA and line of sight 4|$_{.}^{\circ}$|3 4|$_{.}^{\circ}$|9

The GPoA plane parameters using the observed positions and line-of-sight velocities of the M31 satellite galaxies (Conn et al. 2012; McConnachie 2012). They have been determined using the same method applied to the simulated satellite systems. The first column gives the symbol of the quantity described in the second column. The third column assumes that Andromeda XXVII is close to M31, such that it is within the volume covered by the simulated satellite systems. The fourth column assumes that it is outside and should therefore not be considered in the determination of the GPoA parameters that will be compared to the simulated systems.

Table 2.

Parameters of the observed GPoA.

Symbol Description ‘Near’ sample ‘Far’ sample
|$r_\mathrm{M31}^\mathrm{max}$| Maximum allowed satellite radius from M31 300/h kpc 300/h kpc
dAndXXVII Assumed heliocentric distance to Andromeda XXVII ≈800 kpc ≈1250 kpc
Nsat Number of satellites in sample 24 23
NGPoA Number of satellites making up GPoA 15 14
rper rms height of GPoA satellites from best-fitting plane 10.9 kpc 10.9 kpc
rpar rms distance of GPoA satellites projected into best-fitting plane 147.2 kpc 151.0 kpc
rparall rms distance of all Nsat satellites projected into best-fitting plane 136.8 kpc 138.8 kpc
Ncoorb Number of satellites in GPoA which co-orbit 13 13
c/a Minor-to-major axis rms height ratio 0.105 0.102
b/a Intermediate-to-major axis rms height ratio 0.788 0.740
Nlopside Number of GPoA satellites on MW side of M31 13 13
αMW Angle between GPoA and line of sight 4|$_{.}^{\circ}$|3 4|$_{.}^{\circ}$|9
Symbol Description ‘Near’ sample ‘Far’ sample
|$r_\mathrm{M31}^\mathrm{max}$| Maximum allowed satellite radius from M31 300/h kpc 300/h kpc
dAndXXVII Assumed heliocentric distance to Andromeda XXVII ≈800 kpc ≈1250 kpc
Nsat Number of satellites in sample 24 23
NGPoA Number of satellites making up GPoA 15 14
rper rms height of GPoA satellites from best-fitting plane 10.9 kpc 10.9 kpc
rpar rms distance of GPoA satellites projected into best-fitting plane 147.2 kpc 151.0 kpc
rparall rms distance of all Nsat satellites projected into best-fitting plane 136.8 kpc 138.8 kpc
Ncoorb Number of satellites in GPoA which co-orbit 13 13
c/a Minor-to-major axis rms height ratio 0.105 0.102
b/a Intermediate-to-major axis rms height ratio 0.788 0.740
Nlopside Number of GPoA satellites on MW side of M31 13 13
αMW Angle between GPoA and line of sight 4|$_{.}^{\circ}$|3 4|$_{.}^{\circ}$|9

The GPoA plane parameters using the observed positions and line-of-sight velocities of the M31 satellite galaxies (Conn et al. 2012; McConnachie 2012). They have been determined using the same method applied to the simulated satellite systems. The first column gives the symbol of the quantity described in the second column. The third column assumes that Andromeda XXVII is close to M31, such that it is within the volume covered by the simulated satellite systems. The fourth column assumes that it is outside and should therefore not be considered in the determination of the GPoA parameters that will be compared to the simulated systems.

For a consistent comparison, we have applied the same selection criteria to the observed positions of the M31 satellite galaxies. Out of the 27 satellites in the PAndAS sample used in the analyses of Ibata et al. (2013) and Conn et al. (2013), two are outside of the main survey area with 150 kpc radius at M31's distance. In addition, Andromeda XVIII is at a distance of ≈460 kpc from M31 and therefore outside of the volume considered around each of the model host galaxies. It consequently has to be removed from the sample as well. As discussed before, we furthermore discriminate between a near and a far position of Andromeda XXVII, which places it either inside or outside of the considered volume. Therefore we arrive at two sets of parameters, one assuming the distance modulus Richardson et al. (2011) report for Andromeda XXVII (lower limit) and one assuming the distance modulus of Conn et al. (2012) (upper limit).

2.3 Results

The results of our analysis are compiled in Table 3. The table includes only the last three sample prescriptions, the first one is incompatible with the observed satellite distribution and can therefore not be compared with it. The second column of the table lists the criteria according to which a simulated satellite system is considered to be similar to the observed one, which are numbered in the first column. The five criteria are:
\begin{equation} r_{\mathrm{per}}^{\mathrm{sim}} \le r_{\mathrm{per}}^{\mathrm{obs}}, \end{equation}
(1)
\begin{equation} r_{\mathrm{par}}^{\mathrm{sim}} \ge r_{\mathrm{par}}^{\mathrm{obs}}, \end{equation}
(2)
\begin{equation} N_{\mathrm{coorb}}^{\mathrm{sim}} \ge N_{\mathrm{coorb}}^{\mathrm{obs}}, \end{equation}
(3)
\begin{equation} (c/a)^{\mathrm{sim}} \le (c/a)^{\mathrm{obs}}, \end{equation}
(4)
\begin{equation} (b/a)^{\mathrm{sim}} \ge (b/a)^{\mathrm{obs}}. \end{equation}
(5)
The first two criteria are quite similar to the last two because they measure either the absolute or the relative extend of the structure, so it does not make too much sense to combine the four of them. But either of these two combinations (1–2 or 4–5) should be combined with the third criterion, since all three criteria were essential in the discovery of the feature by Ibata et al. (2013). Note again that BB14 have only considered combinations of two of these criteria without providing a motivation for this choice. The combined criteria (numbered 6 to 10 in Table 3) then have to be met simultaneously by a satellite sample. The last three columns in Table 3 give the fractions of systems that contain GPoA-like planes which fulfil the respective criteria. The systems consist of either the 24 brightest satellites from the simulation that lie within the selection volume of which 15 are considered for the plane fit (Psimnear) or the 23 brightest satellites from the simulation that lie within the selection volume of which 14 are considered for the plane fit (Psimfar).
Table 3.

Fractions of simulated satellite systems reproducing the observed GPoA parameters.

# Criterion Psimnear (per cent) Psimfar (per cent) Prandfar (per cent)
1 rpersimrperobs 10.55 15.85 3.80
2 rparsimrparobs 25.58 23.26 21.43
3 NcoorbsimNcoorbobs 4.77 1.79 0.28
4 (c/a)sim ≤ (c/a)obs 10.73 14.03 4.11
5 (b/a)sim ≥ (b/a)obs 17.28 25.59 30.29
6 Criteria 1 & 2 1.49 2.05 0.26
7 Criteria 1 & 3 0.90 0.59 0.01
8 Criteria 4 & 5 0.31 0.77 0.15
9 Criteria 1, 2 & 3 0.17 0.09 <0.007
10 Criteria 3, 4 & 5 0.04 0.04 <0.007
NlopsidesimNlopsideobs 3.77 1.60 0.23
αMWsim ≤ αMWobs 6.83 7.23 6.89
# Criterion Psimnear (per cent) Psimfar (per cent) Prandfar (per cent)
1 rpersimrperobs 10.55 15.85 3.80
2 rparsimrparobs 25.58 23.26 21.43
3 NcoorbsimNcoorbobs 4.77 1.79 0.28
4 (c/a)sim ≤ (c/a)obs 10.73 14.03 4.11
5 (b/a)sim ≥ (b/a)obs 17.28 25.59 30.29
6 Criteria 1 & 2 1.49 2.05 0.26
7 Criteria 1 & 3 0.90 0.59 0.01
8 Criteria 4 & 5 0.31 0.77 0.15
9 Criteria 1, 2 & 3 0.17 0.09 <0.007
10 Criteria 3, 4 & 5 0.04 0.04 <0.007
NlopsidesimNlopsideobs 3.77 1.60 0.23
αMWsim ≤ αMWobs 6.83 7.23 6.89

Fractions of simulated satellite systems fulfilling the criteria numbered in the first and specified in the second column. Psimnear is the fraction of simulated systems fulfilling the criteria for the sample assuming Andromeda XXVII to be close to M31 (compared to the parameters given in column three of Table 2), Psimfar those for the simulated systems selected assuming Andromeda XXVII to be far (compared to the parameters in column four of Table 2). Prandfar is based on the same selection criteria as Psimfar, but uses the randomized satellite positions and velocities. Only the criteria combinations printed in italics can be considered to fully reproduce a GPoA-like structure.

Table 3.

Fractions of simulated satellite systems reproducing the observed GPoA parameters.

# Criterion Psimnear (per cent) Psimfar (per cent) Prandfar (per cent)
1 rpersimrperobs 10.55 15.85 3.80
2 rparsimrparobs 25.58 23.26 21.43
3 NcoorbsimNcoorbobs 4.77 1.79 0.28
4 (c/a)sim ≤ (c/a)obs 10.73 14.03 4.11
5 (b/a)sim ≥ (b/a)obs 17.28 25.59 30.29
6 Criteria 1 & 2 1.49 2.05 0.26
7 Criteria 1 & 3 0.90 0.59 0.01
8 Criteria 4 & 5 0.31 0.77 0.15
9 Criteria 1, 2 & 3 0.17 0.09 <0.007
10 Criteria 3, 4 & 5 0.04 0.04 <0.007
NlopsidesimNlopsideobs 3.77 1.60 0.23
αMWsim ≤ αMWobs 6.83 7.23 6.89
# Criterion Psimnear (per cent) Psimfar (per cent) Prandfar (per cent)
1 rpersimrperobs 10.55 15.85 3.80
2 rparsimrparobs 25.58 23.26 21.43
3 NcoorbsimNcoorbobs 4.77 1.79 0.28
4 (c/a)sim ≤ (c/a)obs 10.73 14.03 4.11
5 (b/a)sim ≥ (b/a)obs 17.28 25.59 30.29
6 Criteria 1 & 2 1.49 2.05 0.26
7 Criteria 1 & 3 0.90 0.59 0.01
8 Criteria 4 & 5 0.31 0.77 0.15
9 Criteria 1, 2 & 3 0.17 0.09 <0.007
10 Criteria 3, 4 & 5 0.04 0.04 <0.007
NlopsidesimNlopsideobs 3.77 1.60 0.23
αMWsim ≤ αMWobs 6.83 7.23 6.89

Fractions of simulated satellite systems fulfilling the criteria numbered in the first and specified in the second column. Psimnear is the fraction of simulated systems fulfilling the criteria for the sample assuming Andromeda XXVII to be close to M31 (compared to the parameters given in column three of Table 2), Psimfar those for the simulated systems selected assuming Andromeda XXVII to be far (compared to the parameters in column four of Table 2). Prandfar is based on the same selection criteria as Psimfar, but uses the randomized satellite positions and velocities. Only the criteria combinations printed in italics can be considered to fully reproduce a GPoA-like structure.

In the following we focus on the ‘far’ sample prescription, which assumes the same distance (≈1250 kpc) to Andromeda XXVII as Ibata et al. (2013) and Conn et al. (2013), but we note that the results for the sample prescription assuming a closer distance to that satellite are similar (see the results compiled in Table 3).

2.3.1 Comparing the plane thickness, radial extend and co-orbiting sense

Fig. 1 plots the rms distance rpar of the satellites associated with the planes versus the rms thickness rper of the fitted planes similar to fig. 2 in BB14. In contrast to these authors we cannot find rpar values much larger than 300 kpc because our method only considers satellites within 426 kpc from the host avoiding satellites with larger distances that are situated in the corners of the cube extracted from the simulation.

Figure 1.

rms distance (rpar) of the satellites projected on to their best-fitting plane plotted against the rms height (rper) of the satellites from the plane for the ‘far’ sample. Satellite systems with a co-orbiting number NcoorbsimNcoorbobs are plotted as red dots, all others as black points. The blue dot denotes the properties of the observed GPoA (determined for the satellite system after being subject to the same selection criteria as the simulated ones, see Table 2 for the parameters). Only satellite systems falling in the area to the left and above the blue dashed lines reproduce the observed situation, and there are only 14 out of 15 000 simulated systems (0.09 per cent, criterion #9 in Table 3) which also contain a sufficient number of co-orbiting satellites. Compare to fig. 2 in BB14.

To be comparable to the observed GPoA, a simulated satellite population must contain a plane of satellites which is similarly narrow. Like BB14 we first measure this using the rms height rper of the fitted satellite plane, i.e. rpersimrperobs. Only 15.9 per cent of the simulated satellite samples have this property. That number is much lower than 40 per cent reported by BB14, even though we are using the same data. The difference indicates that our adjustments to the selection of simulated satellites, i.e. applying the same selection to the observed satellite sample and the simulated satellites, have a significant influence on the result. Only 3.8 per cent of the sample using the randomized satellite positions fulfil this criterion, indicating that the simulated satellite systems are not fully isotropic. This is to be expected, because dark matter haloes and their subhalo systems are not perfectly spherical, due to their merger history and the preferential infall of mass along dark matter filaments (e.g. Libeskind et al. 2011). That some non-isotropic satellite distributions have a larger likelihood to form satellite planes has already been demonstrated by Metz et al. (2007). However, the non-isotropy of the simulated satellite systems is generally insufficient to explain the strong observed coherence of the satellite galaxy positions and velocities (e.g. Pawlowski et al. 2012b).

We find rparsim is larger than rparobs in 23.3 per cent of all simulated satellite samples. This fraction is again significantly different from the 10 per cent reported by BB14. The reason for the discrepancy is that we determine rpar only for those satellites in the GPoA that are close enough to M31 to be within the considered simulation volumes. This either excludes Andromeda XXVII or assumes it to be close to M31. Both scenarios result in a considerably smaller rparobs than the one apparently used by BB14.5 This second criterion is fulfilled by a similar fraction of 21.4 per cent of the samples using randomized satellite positions, as to be expected because these satellite systems have the same radial distribution.

To resemble the planar GPoA, criteria 1 and 2 have to be met simultaneously. Only 2.0 per cent of all simulated satellite systems do so (criterion 6 in Table 3), a value consistent with the 2 per cent reported by BB14. In addition, we know that the GPoA satellites preferentially co-orbit. The number of co-orbiting satellites in the plane, Ncoorb, must be at least as high as observed: NcoorbsimNcoorbobs. This individual criterion (number 3) is met by only 1.8 per cent of the simulated satellite samples. Only 0.28 per cent of the randomized satellite systems fulfil this third criterion.

Only when the third criterion (Ncoorb) is combined with the other two (rper and rpar) can a simulated satellite system be considered similar to the GPoA (although this is still ignoring the alignment of the GPoA with the line of sight and the lopsidedness of the observed satellite distribution). Unfortunately BB14 do not discuss this criteria combination. We find that only 0.09 per cent (14 out of 15 000 samples) meet all three criteria simultaneously (criterion 9 in Table 3).

Finally, none of the 15 000 systems based on the randomized satellite positions reproduce all three criteria. From this we can estimate an upper limit of 0.007 per cent, which is consistent with Ibata et al. (2013), who report a significance of 99.998 per cent for the GPoA, meaning that a fully isotropic satellite distribution can reproduce the observed properties in only 0.002 per cent of all cases.

2.3.2 Comparing the plane axis ratios and co-orbiting sense

The axis ratios of the fitted planes are plotted in Fig. 2. Most simulated satellite systems have a ratio of short-to-long axis (c/a) in the range of 0.1 to 0.3, but they show a wide spread in the intermediate-to-long axis ratio b/a. Most importantly, there is a distinct trend of larger c/a for larger b/a. This indicates that narrow distributions are more prolate while more oblate (and therefore plane-like) distributions are wider along their shortest axis. The danger of misinterpreting elongated, prolate satellite distributions as planes similar to the GPoA when only looking at the rms height is therefore real.

Figure 2.

Axis ratios of the GPoA-like planes in the simulated satellite systems for the ‘far’ sample. Satellite systems with a co-orbiting number NcoorbsimNcoorbobs are plotted as red dots, all others as black points. The blue dot denotes the axis ratios of the observed GPoA (determined for the satellite system after being subject to the same selection criteria as the simulated ones, see Table 2 for the parameters). Only satellite systems falling into the area below and to the right of the blue dashed lines reproduce the observed situation, and there are only six out of 15 000 simulated systems (0.04 per cent, criterion #10 in Table 3) which also contain a sufficient number Ncoorb of co-orbiting satellites. These systems are not the same as the ones fulfilling all three criteria shown in Fig. 1.

For a simulated satellite population to display a similarly narrow alignment as the observed GPoA, (c/a) has to be at least as small as the observed value, i.e. (c/a)sim ≤ (c/a)obs. In addition, for the alignment to be similarly plane-like, b/a has to be at least as large as for the observed system, i.e. (b/a)sim ≥ (b/a)obs. Only 0.8 per cent of the simulated satellite samples fulfil these two criteria simultaneously. In addition, to resemble the co-orbiting satellite plane the Ncoorb criterion has to be fulfilled too. Only six, i.e. 0.04 per cent, of all simulated satellite systems fulfil these three criteria simultaneously and can be considered to resemble the observed GPoA (not accounting for the alignment with the line of sight and the lopsidedness). This number is lower than the one found for the previously discussed three-criteria combination and corresponds to different satellite systems. Again, none of the randomized satellite systems fulfil all three criteria simultaneously.

2.3.3 Effects of sample selection

Fig. 3 compares the cumulative probability distributions of the three quantities rper, rpar and Ncoorb for the four samples as defined in Section 2.2. The major difference is that for the sample prescriptions constructed for consistency with the observed satellite system (near and far samples, ii and iii in Section 2.2), the rms height of the planes is considerably larger than for the sample selection used by BB14. This is because the latter draws a larger number of satellites than observed from a survey volume which is smaller than the observed one. The sample using randomized satellite positions results in even larger rper. Thus, the average rms heights of the fitted planes depend strongly on the sample selection recipe. It is smallest for the BB14-like sample (〈rperBB14 = 12.46 kpc) and larger for the samples containing fewer satellites (〈rpernear = 16.33 and 〈rperfar = 15.15 kpc). The difference between 〈rperBB14 and 〈rperfar is similar to the difference between the latter and 〈rperrand = 17.57. The average axis ratios show the same behaviour.

Figure 3.

Cumulative probability functions of the rms height rper of satellites from the best-fitting plane (left-hand panel), their rms distance rpar when projected on to the plane (central panel) and the number of co-orbiting satellites Ncoorb associated with the plane (right-hand panel). Each panel shows four different distributions illustrating the effects of different sample selection prescription. Dashed blue line: Nsat = 27, NGPoA = 15, artificial PAndAS area is modelled as a cylinder of radius 128 kpc and satellites are drawn from a randomly oriented cube with edge length of 300/h kpc. Solid red line: Nsat = 23, NGPoA = 14. Dotted green line: Nsat = 24, NGPoA = 15. Thin solid grey line: Nsat = 23, NGPoA = 14, but drawing satellites not from the simulation directly but first randomizing their positions and velocities. In the latter three cases the artificial PAndAS survey is modelled as a cone with a radius of 150 kpc at the distance of M31 and the satellites are drawn from a randomly oriented sphere with radius 300/h kpc. It is more likely to find narrow planes (small rper) when more satellite galaxies are drawn from a more narrow volume. Drawing sample galaxies from a cube produces more satellite systems with large rpar values, as sometimes the random line of sight is oriented along two opposing corners. In satellite planes consisting of 14 instead of 15 objects the probability of finding large numbers of co-orbiting satellites is lower.

In addition, for the sample prescriptions constructed to be consistent with the observed satellite sample, the probability of finding large values of rpar is lower while that of finding small values is increased. This clearly shows the effect of choosing satellites either from the randomly oriented cube or from a sphere embedded in them. However, the average rms distance of the GPoA satellites projected into the plane does not depend strongly on the sample selection, the averages are similar for all samples: 〈rparBB14 = 125 kpc, 〈rparnear = 129 kpc, 〈rparfar = 129 kpc and 〈rparrand = 127 kpc.

As expected the number of co-orbiting satellites tends to be lower for the sample considering only NGPoA = 14 satellites to belong to the GPoA instead of 15. It is also much lower for the sample using randomized satellite positions and velocity directions. The average number of co-orbiting satellites in the GPoA-like planes depends mostly on the number of satellites in the GPoA. The samples using 15 GPoA members result have 〈NcoorbBB14 = 9.56 and 〈Ncoorbnear = 9.57, while the one using only 14 satellites for the plane fit has 〈Ncoorbfar = 8.95. Randomizing the satellite positions and velocity directions reduces this number by about 0.5, to 〈Ncoorbrand = 8.48.

2.3.4 Additional properties

In addition to the previously discussed parameters, the GPoA exhibits two more puzzling properties. It is seen almost edge-on (the angle αMW between the line of sight and the GPoA is small) and the distribution of satellites belonging to the GPoA is lopsided, i.e. a significantly larger number (13 of the 15 in the GPoA) of satellite galaxies (Nlopside) is on the near side of M31. Among the simulated satellite samples, 1.6 per cent show a similar degree of lopsidedness (NlopsidesimNlopsideobs) and 7.2 per cent align similarly well with the MW (αMWsim ≤ αMWobs). In combination, these two criteria are met by 0.13 per cent of the samples and combining them with either of the previous three-criteria combinations results in no sample fulfilling the five criteria simultaneously. However, we note that these two features might be related to the dynamics of the LG and the MW's particular viewing position which are not considered in the selection of host galaxies and their viewing directions.

2.4 Summary

A satellite galaxy configuration as observed in the GPoA around M31, satisfying simultaneously the three criteria defined by BB14, is only found in 0.04 to 0.17 per cent of the satellite samples based on the Millennium II simulation, depending on Andromeda XXVII's distance to M31. Additionally to the more realistic sample selection, the main reason for this difference of our results from those of BB14 is that we require GPoA-like planes to agree with all three observed properties, instead of arbitrarily combining only two of them.

Our result is robust against changing the prescription of the major GPoA properties. Comparing the rms height rper, rms distance of the satellites projected into the plane rpar and their co-orbiting number Ncoorb yields similarly low probabilities as a comparison of the axis ratios c/a and b/a in combination with Ncoorb.

Further improvements of the analysis are expected from (i) using higher resolution simulations so the analysis does not have to rely on orphan galaxies and (ii) a comparison of the observed satellite distribution of large galaxies in LG-like galaxy groups. These points will be addressed in a future paper.

While our work was being refereed, Ibata et al. (2014) have also published a re-analysis of the work by BB14. While our analysis stays as close as possible to the one performed by BB14, the analysis of Ibata et al. (2014) remains as close as possible to the one performed on the observational data by Ibata et al. (2013). They apply the exact PAndAS survey area to the simulations, enforce the satellite planes to be seen edge-on as is the case for observed M31 system and in addition check whether the specific angular momenta of the simulated satellites in the GPoA-like planes are as large as observed. It is reassuring that the results of both our studies with different approaches agree: Ibata et al. (2014) find that only 0.04 per cent of the planes in satellite systems drawn from the Millennium II simulation can reproduce the thickness, radial extend and velocity coherence of the GPoA, to be compared with our result of 0.04 to 0.17 per cent.

3 SEARCHING FOR A MW VPOS-LIKE STRUCTURE IN ΛCDM SIMULATIONS

3.1 Previous method and its limitations

W13 report that 6–13 per cent of satellite systems in cosmological simulations are as flat as that of the MW VPOS. However, there are a number of indications that this overestimates the probability of VPOS-like structures to appear in ΛCDM simulations. The major points and our adjustments to the analysis are given below.

  • Arguably the most important property of the observed MW satellite system is that a significant fraction of the satellites co-orbit within the VPOS. This orbital alignment indicates that it is rotationally stabilized and therefore could be long lived (Lynden-Bell 1976; Metz et al. 2008; Pawlowski et al. 2012a; Pawlowski & Kroupa 2013). W13 have not tested whether the simulated satellite systems show a similar degree of phase-space correlation. Because the velocities of satellite galaxies in simulations are known, we have to closely pay attention to this observational constraint on the VPOS properties and will therefore incorporate it in the following analysis.

  • W13 analyse satellite systems in the large-scale Millennium II simulation as well as those around the six main haloes in the high-resolution Aquarius simulations (Springel et al. 2008). None of the Aquarius satellite systems they analyse can reproduce the same amount of flattening as measured for the VPOS using either c/a (their fig. 7) or rper (their fig. 10). As noted by W13, their conclusion that 6–13 per cent of the satellite systems can be as flattened as the VPOS is only based on satellite systems constructed by applying a semi-analytic galaxy formation model to the Millennium II simulation (Guo et al. 2011). This might indicate that the inclusion of orphan galaxies with incorrect positions affects the analysis, but the reason could also be small number statistics as there are only six host galaxies in the Aquarius simulation suite. In addition, W13 use satellites obtained from the dark-matter-only simulations using slightly different semi-analytical models for the Aquarius (Cooper et al. 2010) and Millennium II (Guo et al. 2011) simulations, while in addition using an abundance matching technique for the Aquarius simulations, too. In the following we will focus on the Millennium II simulation.

  • W13 compare the observed and simulated satellite systems using two different measures of the flattening: the axis ratios of the satellite distribution c/a, where the observed value is met by 6 per cent of their simulated systems, and the ratio of the rms thickness rper of the best-fitting plane to a fixed cut-off radius of rcut = 250 kpc. As already pointed out in Pawlowski et al. (2013), the latter method, according to which 13 per cent of the simulated systems are meeting the observed value, is biased. It does not consider the radial distribution of the satellites, and thus a radially more concentrated simulated satellite system could count as VPOS-like even though the satellites are distributed nearly isotropically. We will account for this by using the same method as in Section 2 and demanding the satellites to be sufficiently extended radially, as measured by rpar.

  • Another important aspect influencing the result of W13 is how they incorporate the obscuration by the Galactic disc in their satellite selection. They assume that no satellite galaxy can be observed within a critical angular distance from the MW disc (θcrit = 9|$_{.}^{\circ}$|5), and model this selection effect by placing obscuring discs into their simulated satellite systems, ignoring all satellites within θcrit of them. Instead of using randomly oriented obscuring discs they iteratively orient them such that they are perpendicular to the major axis of the satellite distribution. While the observed VPOS is indeed almost perpendicular to the MW plane, there is no reason to expect such an orientation a priori. In fact, observations and cosmological simulations suggest the contrary, i.e. a preferred planar alignment of satellites with the major axis of their host galaxies (e.g. Sales & Lambas 2004, 2009; Brainerd 2005; Bett et al. 2010; Lovell et al. 2011). Iteratively forcing the obscuring disc to be perpendicular to the major axis of the satellite distribution has the effect of obscuring those satellite galaxies which are perpendicular (and therefore on average most distant) from the best-fitting satellite plane. This artificially increases the measured flattening of the system, biasing them to find an agreement with the observed, narrow VPOS. If the reason for the VPOS is the absence of satellites in the disc of the MW (Zone of Avoidance), then incorporating randomly oriented obscuring discs already accounts for the effect of the MW disc on the satellite distribution and should result in finding plenty of VPOS-like structures in the simulated satellite systems.

  • When determining the properties of the observed MW satellite system which will be compared with those of the simulated systems, W13 decide to use the 11 brightest, i.e. classical MW satellite galaxies. They note that the 12th-brightest MW satellite Canes Venatici I (CVnI) has only a slightly lower luminosity than the next-brighter classical satellites Draco and Ursa Minor (MV = −8.6 compared to −8.8; McConnachie 2012), but argue that the inclusion of CVnI in the MW satellite sample has almost no effect on the axis ratio of the satellite system and do not include it in their analysis. However, one of their results is that simulated satellite systems become less flattened when more satellites are considered (see e.g. their fig. 9). Thus, even though the flattening might be the same for 11 or 12 observed satellites, there might be a difference between the flattening of 11 and 12 simulated satellites. Whether the addition of one object changes the results might be doubted, but should be tested. We will do so in the following.

Following the same procedure as in the previous section, we account for these issues with a refined analysis method and compare the results to those which W13 obtained using their method.

3.2 New method

For the following analysis we make use of the same data from the Millennium II simulation as in Section 2. This results in a similar main halo selection as used by W13, but in a slightly smaller range in virial mass of 1.1 to 1.7 × 1012 M instead of 1.0 to 2.0 × 1012 M. Working with the same main halo sample has the advantage of a direct comparison of the probabilities for finding either a VPOS or a GPoA-like system among the same main galaxy population. This makes it possible to decide which of the structures is the least common in the ΛCDM satellite systems we analyse. By combining the probabilities of finding a VPOS or a GPOA we are also able to determine the probability to find both structures in the LG, under the reasonable assumption that both the MW and M31 have halo masses within the same range.

We only consider satellites with distances from their host in the range of 15 to 260 kpc. The upper limit was chosen to be compatible with the Galactocentric distance of Leo I, the most distant of the 11 classical satellite galaxies. For each main halo, we select the 11 satellite galaxies with the largest baryonic masses and with an angular separation θ > θcrit from a randomly oriented plane, imitating the obscuration of the Galactic disc. We chose the same obscuration angle as W13, θcrit = 9|$_{.}^{\circ}$|5. However, in contrast to their analysis we do not require that this galactic plane is perpendicular to the minor axis of the plane fitted to the satellite galaxies, for reasons we discussed before.

In total there are 1693 main haloes in the chosen mass range which contain a sufficient number of satellite galaxies. As in Section 2 we generate 10 samples for each main halo, by adopting different random obscuring disc orientations, such that we have 16 930 simulation samples.

For each sample we determine the best-fitting plane of all satellites using the method described in Section 2. Again we compute the relevant properties (in particular rper, rpar, c/a  and  b/a) and compare them to those derived for the observed MW satellite galaxy system. Further, we determine how strongly the satellites are co-orbiting employing the method presented in Pawlowski et al. (2012b). For this purpose, we construct the orbital poles |$\hat{n_i}$| (directions of angular momenta) of the simulated satellite galaxies and determine which of all possible combinations of eight of these poles for a given satellite system have the smallest spherical standard distance Δstd (Metz et al. 2007; Pawlowski et al. 2012b) from their average direction |$\langle \hat{n} \rangle$|⁠,
\begin{equation*} \Delta _{\mathrm{std}} = \sqrt{ \frac{ \sum _{i=1}^{8} \left[ \arccos \left( \langle \hat{n} \rangle \cdot \hat{n}_{i} \right) \right] ^{2} }{ 8 } } = \min . \end{equation*}
These eight orbital poles are the best-aligned (or best-concentrated) orbital poles of the sample, and we record their concentration as measured by the spherical standard distance Δstd as well as the angle between their average direction and the normal vector of the best-fitting plane θVPOS. As in Pawlowski & Kroupa (2013) we checked that considering only the six best-aligned poles does not change the result.

In addition to the Millennium II data we also compare our results to randomized satellite systems. For these, we rotate both the positions and velocities of the satellite galaxies into random directions before applying the obscuration cut. This process leads to isotropic satellite systems with the same radial distribution as in the non-randomized case.

To test the statistical robustness of our results we also select systems of 12 satellite galaxies to be compared to the properties of plane fitted to the 11 classical MW satellite galaxies plus CVnI, which has a very similar luminosity as the least-luminous classical satellites Draco and Ursa Minor (as discussed in 3.1). Because CVnI was discovered in the Sloan Digital Sky Survey (SDSS) survey at large Galactic latitude (b = +79|$_{.}^{\circ}$|8), we require the 12th satellite galaxy drawn from the simulations to be within 60° from the north pole direction of the artificial obscuring disc.

Table 4 compiles the properties of the VPOS planes fitted to the observed MW satellite galaxies for samples consisting of the 11 classical satellite galaxies and for the sample extended with CVnI. The orbital pole of CVnI is not known, so we adopt the same orbital pole concentration and average direction as determined from the smaller sample of only 11 satellites. We will compare these properties with those of the larger simulation samples, i.e. we will compare the concentration of the eight best-aligned orbital poles out of only 11 observed satellite galaxies with the concentration of the eight best-aligned orbital poles out of 12 poles in the simulation. The benefit is that our results will remain valid once CVnI's proper motion is known, even in the unlikely case that it is not co-orbiting in the VPOS. Our results are therefore lower limits, because assuming it to orbit in VPOS, CVnI's predicted proper motion would make it align better than some of the known poles (Pawlowski & Kroupa 2013).

Table 4.

Parameters of the observed VPOS.

Quantity Description Classical sats. Classical sats. + CVnI
NVPOS Number of satellites used for the VPOS plane fit 11 12
rper rms height of satellites from best-fitting plane 19.6 kpc 20.8 kpc
rpar rms distance of satellites projected into best-fitting plane 129.5 kpc 138.9 kpc
c/a Minor-to-major axis rms height ratio 0.182 0.183
b/a Intermediate-to-major axis rms height ratio 0.508 0.530
Δstd Spherical standard distance of eight most concentrated orbital poles 29|$_{.}^{\circ}$|3 ≤29|$_{.}^{\circ}$|3
θVPOS Angle between VPOS normal and average direction of eight most concentrated orbital poles 18|$_{.}^{\circ}$|9 24|$_{.}^{\circ}$|5
Quantity Description Classical sats. Classical sats. + CVnI
NVPOS Number of satellites used for the VPOS plane fit 11 12
rper rms height of satellites from best-fitting plane 19.6 kpc 20.8 kpc
rpar rms distance of satellites projected into best-fitting plane 129.5 kpc 138.9 kpc
c/a Minor-to-major axis rms height ratio 0.182 0.183
b/a Intermediate-to-major axis rms height ratio 0.508 0.530
Δstd Spherical standard distance of eight most concentrated orbital poles 29|$_{.}^{\circ}$|3 ≤29|$_{.}^{\circ}$|3
θVPOS Angle between VPOS normal and average direction of eight most concentrated orbital poles 18|$_{.}^{\circ}$|9 24|$_{.}^{\circ}$|5

The VPOS plane parameters using the observed positions and line-of-sight velocities of the brightest MW satellite galaxies using the data compiled by McConnachie (2012). They have been determined using the same method applied to the simulated satellite systems. The first column gives the symbol of the quantity described in the second column. The third column gives the parameters of the plane fitted to the 11 classical satellite galaxies, while the fourth column gives those of the plane fitted to these satellites plus the next-fainter satellite CVnI.

Table 4.

Parameters of the observed VPOS.

Quantity Description Classical sats. Classical sats. + CVnI
NVPOS Number of satellites used for the VPOS plane fit 11 12
rper rms height of satellites from best-fitting plane 19.6 kpc 20.8 kpc
rpar rms distance of satellites projected into best-fitting plane 129.5 kpc 138.9 kpc
c/a Minor-to-major axis rms height ratio 0.182 0.183
b/a Intermediate-to-major axis rms height ratio 0.508 0.530
Δstd Spherical standard distance of eight most concentrated orbital poles 29|$_{.}^{\circ}$|3 ≤29|$_{.}^{\circ}$|3
θVPOS Angle between VPOS normal and average direction of eight most concentrated orbital poles 18|$_{.}^{\circ}$|9 24|$_{.}^{\circ}$|5
Quantity Description Classical sats. Classical sats. + CVnI
NVPOS Number of satellites used for the VPOS plane fit 11 12
rper rms height of satellites from best-fitting plane 19.6 kpc 20.8 kpc
rpar rms distance of satellites projected into best-fitting plane 129.5 kpc 138.9 kpc
c/a Minor-to-major axis rms height ratio 0.182 0.183
b/a Intermediate-to-major axis rms height ratio 0.508 0.530
Δstd Spherical standard distance of eight most concentrated orbital poles 29|$_{.}^{\circ}$|3 ≤29|$_{.}^{\circ}$|3
θVPOS Angle between VPOS normal and average direction of eight most concentrated orbital poles 18|$_{.}^{\circ}$|9 24|$_{.}^{\circ}$|5

The VPOS plane parameters using the observed positions and line-of-sight velocities of the brightest MW satellite galaxies using the data compiled by McConnachie (2012). They have been determined using the same method applied to the simulated satellite systems. The first column gives the symbol of the quantity described in the second column. The third column gives the parameters of the plane fitted to the 11 classical satellite galaxies, while the fourth column gives those of the plane fitted to these satellites plus the next-fainter satellite CVnI.

3.3 Results

The results of our analysis are compiled in Table 5. As in Section 2, the second column of the table lists the criteria according to which a simulated satellite system is considered to be similar to the observed one, numbered in the first column. Combined criteria (numbers 7 to 17) have to be met simultaneously. The last three columns give the fractions of systems fulfilling the respective criteria which consist of either the 11 brightest satellites in the simulation outside of the obscuration region (Psim11), the 11 brightest outside the obscuration area for randomized angular positions and velocities of the satellites (Prand11) and for the 12 brightest satellites in the simulation and outside the obscured region, of which the least-luminous one lies within the approximate SDSS footprint (Psim12).

Table 5.

Fractions of simulated satellite systems reproducing the observed VPOS parameters.

# Criterion Psim11 (per cent) Prand11 (per cent) Psim12 (per cent)
1 rpersimrperobs 4.45 2.10 4.44
2 rparsimrparobs 51.65 49.69 34.90
3 (c/a)sim ≤ (c/a)obs 3.35 0.92 2.39
4 (b/a)sim ≥ (b/a)obs 86.36 92.86 84.05
5 Δstdsim ≤ Δstdobs 2.48 0.16 4.11
6 θVPOSsim ≤ θVPOSobs 13.44 10.06 24.11
7 Criteria 1 & 2 1.22 0.27 0.63
8 Criteria 3 & 4 1.57 0.64 1.02
9 Criteria 1 & 5 0.35 <0.006 0.41
10 Criteria 3 & 5 0.28 0.01 0.27
11 Criteria 5 & 6 0.97 0.05 1.92
12 Criteria 1, 2 & 5 0.15 <0.006 0.11
13 Criteria 3, 4 & 5 0.12 <0.006 0.08
14 Criteria 1, 5 & 6 0.15 <0.006 0.23
15 Criteria 3, 5 & 6 0.12 0.01 0.15
16 Criteria 1, 2, 5 & 6 0.06 <0.006 0.04
17 Criteria 3, 4, 5 & 6 0.05 <0.006 0.05
# Criterion Psim11 (per cent) Prand11 (per cent) Psim12 (per cent)
1 rpersimrperobs 4.45 2.10 4.44
2 rparsimrparobs 51.65 49.69 34.90
3 (c/a)sim ≤ (c/a)obs 3.35 0.92 2.39
4 (b/a)sim ≥ (b/a)obs 86.36 92.86 84.05
5 Δstdsim ≤ Δstdobs 2.48 0.16 4.11
6 θVPOSsim ≤ θVPOSobs 13.44 10.06 24.11
7 Criteria 1 & 2 1.22 0.27 0.63
8 Criteria 3 & 4 1.57 0.64 1.02
9 Criteria 1 & 5 0.35 <0.006 0.41
10 Criteria 3 & 5 0.28 0.01 0.27
11 Criteria 5 & 6 0.97 0.05 1.92
12 Criteria 1, 2 & 5 0.15 <0.006 0.11
13 Criteria 3, 4 & 5 0.12 <0.006 0.08
14 Criteria 1, 5 & 6 0.15 <0.006 0.23
15 Criteria 3, 5 & 6 0.12 0.01 0.15
16 Criteria 1, 2, 5 & 6 0.06 <0.006 0.04
17 Criteria 3, 4, 5 & 6 0.05 <0.006 0.05

Fractions of simulated satellite systems fulfilling the criteria numbered in the first column and specified in the second. Psim11 gives the fraction of simulated systems fulfilling the criteria for the 11 brightest satellites not hidden by the obscuring MW disc, Prand11 those for the 11 brightest satellites with randomized position and velocity directions (both compared to the parameters given in column three of Table 4). Psim12 gives the fractions for the simulated satellite systems consisting of the 11 brightest satellites outside the obscured region and the next-fainter one within the approximate SDSS survey area (i.e. compared to the parameters determined including CVnI given in column four of Table 4). Only the criteria combinations printed in italics can be considered to fully reproduce a VPOS-like structure.

Table 5.

Fractions of simulated satellite systems reproducing the observed VPOS parameters.

# Criterion Psim11 (per cent) Prand11 (per cent) Psim12 (per cent)
1 rpersimrperobs 4.45 2.10 4.44
2 rparsimrparobs 51.65 49.69 34.90
3 (c/a)sim ≤ (c/a)obs 3.35 0.92 2.39
4 (b/a)sim ≥ (b/a)obs 86.36 92.86 84.05
5 Δstdsim ≤ Δstdobs 2.48 0.16 4.11
6 θVPOSsim ≤ θVPOSobs 13.44 10.06 24.11
7 Criteria 1 & 2 1.22 0.27 0.63
8 Criteria 3 & 4 1.57 0.64 1.02
9 Criteria 1 & 5 0.35 <0.006 0.41
10 Criteria 3 & 5 0.28 0.01 0.27
11 Criteria 5 & 6 0.97 0.05 1.92
12 Criteria 1, 2 & 5 0.15 <0.006 0.11
13 Criteria 3, 4 & 5 0.12 <0.006 0.08
14 Criteria 1, 5 & 6 0.15 <0.006 0.23
15 Criteria 3, 5 & 6 0.12 0.01 0.15
16 Criteria 1, 2, 5 & 6 0.06 <0.006 0.04
17 Criteria 3, 4, 5 & 6 0.05 <0.006 0.05
# Criterion Psim11 (per cent) Prand11 (per cent) Psim12 (per cent)
1 rpersimrperobs 4.45 2.10 4.44
2 rparsimrparobs 51.65 49.69 34.90
3 (c/a)sim ≤ (c/a)obs 3.35 0.92 2.39
4 (b/a)sim ≥ (b/a)obs 86.36 92.86 84.05
5 Δstdsim ≤ Δstdobs 2.48 0.16 4.11
6 θVPOSsim ≤ θVPOSobs 13.44 10.06 24.11
7 Criteria 1 & 2 1.22 0.27 0.63
8 Criteria 3 & 4 1.57 0.64 1.02
9 Criteria 1 & 5 0.35 <0.006 0.41
10 Criteria 3 & 5 0.28 0.01 0.27
11 Criteria 5 & 6 0.97 0.05 1.92
12 Criteria 1, 2 & 5 0.15 <0.006 0.11
13 Criteria 3, 4 & 5 0.12 <0.006 0.08
14 Criteria 1, 5 & 6 0.15 <0.006 0.23
15 Criteria 3, 5 & 6 0.12 0.01 0.15
16 Criteria 1, 2, 5 & 6 0.06 <0.006 0.04
17 Criteria 3, 4, 5 & 6 0.05 <0.006 0.05

Fractions of simulated satellite systems fulfilling the criteria numbered in the first column and specified in the second. Psim11 gives the fraction of simulated systems fulfilling the criteria for the 11 brightest satellites not hidden by the obscuring MW disc, Prand11 those for the 11 brightest satellites with randomized position and velocity directions (both compared to the parameters given in column three of Table 4). Psim12 gives the fractions for the simulated satellite systems consisting of the 11 brightest satellites outside the obscured region and the next-fainter one within the approximate SDSS survey area (i.e. compared to the parameters determined including CVnI given in column four of Table 4). Only the criteria combinations printed in italics can be considered to fully reproduce a VPOS-like structure.

3.3.1 A thin VPOS

W13 have reported that 6 per cent of the simulated satellite systems have a smaller axis ratio c/a than observed (our criterion 3 in Table 5). For the observed MW satellites, they measure a value of c/a = 0.18, while we use c/a = 0.182. We therefore use the same criterion, but find only a fraction of 3.4 per cent of the samples to fulfil this criterion. In addition to this relative shape of the distribution as measured by the axis ratio, W13 also test how often their satellite systems are as thin as the observed one by using the rms height rper of the satellites from the best-fitting plane as measured in kpc. They compare the ratio of rper to a fixed cut-off radius of Rcut = 250 kpc, chosen to be similar to the distance of the most distant known MW satellite Leo I. As Rcut is kept at a fixed value independent of the simulated satellite systems, this method is identical to comparing the rper values of the observed and the simulated satellite systems directly, corresponding to our criterion 1 in Table 5. W13 give a value of |$\frac{r_{\mathrm{per}}}{R_{\mathrm{cut}}} = 0.074$|⁠, which corresponds to rper = 18.5 kpc, whereas we find rper = 19.6 kpc. Our criterion is therefore slightly less strict. In contrast to W13, who report 13 per cent of their systems are having a thickness similar to the observed one, we find this situation is encountered in only 4.5 per cent of our samples. Part of these discrepancies will be due to the enforced polar alignment of their satellite planes with respect to their obscuring disc, which makes sure the part that is most distant from the best-fitting plane is effectively hidden in the obscuration volume. The other obvious difference between the samples is that we use the Millennium II simulation scaled to cosmological parameters from WMAP7 (instead of WMAP1), while the semi-analytic modelling is almost identical (Guo et al. 2013). Thus, using only one of the single-criterion comparisons (c/a or rper) it can already be concluded now that it is considerably less likely for a simulated satellite system to reproduce the observed flattening VPOS than claimed by W13, if the MW discs are oriented randomly.

As discussed in Pawlowski et al. (2013), a more concentrated satellite distribution will naturally have a smaller rms height rper on average. This can be seen in Fig. 4, which plots rper against rpar, the rms distance of satellites projected on to the best-fitting plane: small rper are preferentially present in satellite systems with small rpar. Therefore, following the criteria of BB14, in order to resemble the spatial distribution of the MW satellites the simulated systems not only have to be narrow (rpersimrperobs) but should also be at least as radially extended (rparsimrparobs). This is the case for only 1.2 per cent of the systems.

Figure 4.

rms distance rpar of the satellites projected into their best-fitting plane plotted against the rms height rper of the satellites from the plane for the samples of 11 satellites. Those satellite systems which have a sufficiently concentrated orbital pole distribution Δstdsim ≤ Δstdobs are plotted as red dots, all others as black points. The blue dot denotes the axis ratios of the observed classical MW satellites. Only satellite systems falling into the area to the left and above of the blue dashed lines reproduce the observed spatial distribution. Only 1.22 per cent of the sample systems fall into this range (criterion #7 in Table 5), and only 0.15 per cent simultaneously have sufficiently concentrated orbital poles (criterion #12 in Table 5).

Requiring the relative shape of the satellite distribution to be at least as flattened as observed ((c/a)sim ≤ (c/a)obs and (b/a)sim ≥ (b/a)obs) results in a similar fraction of 1.6 per cent of the systems (see Fig. 5). If the satellite positions are made isotropic by randomizing their angular positions, the probability to find a system similar to the observed one is lower (0.3 and 0.6 per cent for the combined criteria, respectively).

Figure 5.

Axis ratios of the planes fitted to the simulated satellite samples for the samples of 11 satellites. Those satellite systems which have a sufficiently concentrated orbital pole distribution Δstdsim ≤ Δstdobs are plotted as red dots, all others as black points. The blue dot denotes the axis ratios of the observed classical MW satellites. Only satellite systems falling into the area to the right and below by the blue dashed lines reproduce the observed spatial flattening. Only 1.57 per cent of the sample system fall into this range (criterion #8 in Table 5) and only 0.12 per cent simultaneously have sufficiently concentrated orbital poles (criterion #13 in Table 5). Compare to fig. 8 in W13.

3.3.2 A co-orbiting VPOS

Aside from the spatial distribution, a major feature of the VPOS is the co-orbiting sense of the satellite galaxies. This can be expressed by two parameters: the spherical standard distance of the eight best-aligned orbital poles, Δstd, and the alignment of their average direction with the normal vector to the best-fitting satellite plane, θVPOS. Applying this test to the simulated satellite systems results in criteria 5 and 6 in Table 5. In 2.5 per cent of the simulated satellite systems eight of the 11 satellites have orbital poles that align to better than Δstdobs = 29|$_{.}^{\circ}$|3 determined from observational data. This is the case for only 0.16 per cent in the randomized satellite sample, indicating that there is indeed some degree of coherence in the orbital direction of dark matter subhaloes. These fractions are comparable to those reported in Pawlowski & Kroupa (2013). For a random distribution of orbital poles they find a probability of 0.10 per cent that Δstdsim ≤ Δstdobs, while for subhaloes drawn at random from the Aquarius or Via Lactea simulations they find an average of 0.64 per cent and a maximum value of 1.67 per cent for the Aquarius halo C2. The analysis of Pawlowski & Kroupa (2013) does not exclude objects close to the MW plane, which explains why the fractions we find now are somewhat larger. The average direction of the best-aligned orbital poles aligns with the normal direction of the best-fitting plane in 13.4 per cent of the simulated samples and in 10.1 per cent of the randomized ones. The agreement of these values indicates that there is only a weak tendency for the simulated galaxies to align with the plane defined by their positions.

To reproduce the situation of satellite galaxies co-orbiting in the VPOS, the simulated satellite orbital poles have to be both concentrated and close to the normal defined by the plane fitted to their positions, i.e. they have to fulfil both criteria simultaneously. Only 1.0 per cent of the simulated satellite systems do so while such a situation is essentially non-existent for the randomized satellite distributions (upper limit 0.006 per cent).

3.3.3 A thin and co-orbiting VPOS

When considering the spatial distribution or the orbital pole distribution separately only about 1 per cent of the simulated satellite systems can reproduce the observed properties of the MW satellite system. However, to be in agreement with the VPOS both the position as well as the orbital poles have to be considered simultaneously. For both combinations of four criteria (either rper, rpar, Δstd and θVPOS or c/a, b/a, Δstd and θVPOS), only 9, respectively, 8 out of 16 930 simulated satellite systems reproduce the observed situation, only about 0.05 per cent. This is equivalent to less than one of the 1693 host haloes. Even requiring that only three of the criteria are met simultaneously (ignoring either the requirement for the radial distribution of the satellites, rpar, or the alignment of the average orbital pole with the plane normal, θVPOS) does not increase these numbers substantially (criteria 12 to 15 in Table 5). It is in fact already very unlikely that a simulated satellite system happens to have a sufficiently concentrated distribution of orbital poles and a small rms height (0.35 per cent, criterion 9, see Fig. 6) or axis ratio (0.28 per cent, criterion 10) at the same time.

Figure 6.

Orbital pole concentration Δstd plotted against rper, the rms height of the satellite galaxies from the plane fitted to the positions of the 11 brightest satellites around each host. Shown are the results for the samples of 11 satellites. The black points represent the results for the Millennium II simulations, the blue dot the values for the observed satellites around the MW. To be as narrow and as coherently rotating as the MW satellite system, the simulations must lie in the region to the left and below the dashed lines. This is the case for only 0.35 per cent of all simulated satellite samples (criterion #9 in Table 5). If the average orbital pole is simultaneously required to be as close to the plane normal as is the case for the observed VPOS, only 0.15 per cent of the simulated satellite systems can reproduce the observed situation (criterion #14 in Table 5).

3.3.4 Including Canes Venatici I in the VPOS

W13 arbitrarily excluded the MW satellite CVnI from their sample of observed galaxies. However, they also showed that as the number of simulated satellites in a system increases, the axis ratio c/a increases as well. Our results indicate that this trend has already an effect if only one additional satellite galaxy is included. The fraction of simulated systems fulfilling criterion 3 in Table 5 changes from Psim11 = 3.4 per cent to Psim12 = 2.4 per cent. As discussed above, it is more likely to draw 8 closely aligned orbital poles out of a sample of size 12 instead of size 11. This is reflected in the slightly higher probabilities Psim12 compared to Psim11 to fulfil criterion 5 in Table 5. This probability and those based on it are upper limits and will drop if the orbital pole of CVnI turns out to be well aligned with the others. In any case, the probabilities to reproduce the properties of the observed MW system are similar and close to zero for both sample sizes. Including CVnI in the analysis therefore does not change our main result and slightly reduces the likelihood for a simulated satellite system to resemble the observed VPOS.

3.4 Summary

Whether considering the 11 or 12 MW satellites with the largest baryonic mass, in both cases the observed thinness or flattening, orbital pole coherence and alignment of orbits with the best-fitting plane is reproduced by ≲0.2 per cent of the simulated satellite systems, and less than 0.06 per cent of all satellite systems have the four observed properties simultaneously. The observed MW satellite galaxies constituting the VPOS are simply too well aligned in phase space to be reproduced by subhalo based satellite galaxies found in the Millennium II simulation.

If additional satellite objects were to be considered, such as the faint satellite galaxies (Kroupa et al. 2010), globular clusters and stellar and gaseous streams in the MW halo (Pawlowski et al. 2012a), the probability to find the observed situation in these ΛCDM simulations would be even smaller. It is worth noticing that this analysis relies only on the assumption that all classical MW satellites, brighter than MV ≈ −8 and not obscured by our position within the stellar disc, have been discovered within a radius of 260 kpc. It does not rely on the distribution of rogue ultrafaint dwarf spheroidal galaxies that might be found outside the VPOS in the future.

The results of this section have been obtained by following the method of W13, except for the changes listed in Section 3.2. It would of course be interesting to know if these results can be reproduced if instead of the haloes with largest mass in baryons Mbaryon those with the highest V-band luminosity LV or the highest dark matter mass Mpeak or the largest circular velocity Vpeak before infall would be selected, but this is beyond the scope of this paper and left for future studies. One would expect that the haloes with the largest Mbaryon have a high chance also to be among the ones with the highest LV, Mpeak, Vpeak or other parameters of a similar type. This expectation is supported by the comparison for subhaloes from the Aquarius simulation presented in W13. They find that a fraction of about 80 to 85 per cent of the ≈11 highest ranking satellites are identical to those selected using semi-analytic modelling if the selection parameter is changed to Vpeak or Mpeak. Thus, it appears likely that similar results will be obtained when using these different parameters.

Again it would be helpful to have higher resolution simulations available such that no orphan galaxies have to be included in the analysis. Such studies are currently carried out.

4 THE GPOA VIA COLD STREAM ACCRETION?

4.1 Previous method and its limitations

GB136 claimed that the GPoA can be understood if all satellite galaxies are accreted along cold streams that have been proposed as the dominant accretion mechanism for baryons at high redshift (z ≥ 2) based on theoretical work and simulations (e.g. Dekel et al. 2009). GB13 assume that the satellite galaxies are accreted in narrow streams of cold gas, which they model as straight lines crossing the centre of the host galaxy. To determine how likely it is that a certain fraction of satellites lie within a narrow plane, they then determine how often a fraction of lines close to the fraction of M31 satellites in the GPoA (15 out of 27) of such randomly oriented lines are found in a common plane with a sufficiently small opening angle. Their main result is that it is very likely that two to four out of three to eight randomly oriented lines align within a narrow plane.

There are a number of limitations to the analysis of GB13 and its accord with the observed situation. The most important issues are as follows.

  • The validity and importance of the cold stream accretion mode itself as well as its connection to satellite accretion has not yet been observationally confirmed. In particular, it is known that dwarf spheroidal galaxies have higher than Newtonian internal velocity dispersions, which in the standard context means that they have to be embedded in dark matter haloes. How these subhaloes share exactly the dynamics of cold gas along these thin streams in order to then accrete and form satellite galaxies within dark matter haloes still has to be demonstrated. In addition, recent numerical simulations suggest that the cold flow accretion mode might have been overestimated and could therefore to a certain extend be numerical artefacts of the simulation codes used (Nelson et al. 2013).

  • Cold stream accretion is believed to be most important at high redshifts (z ≳ 2) (e.g. Dekel & Birnboim 2006; Ocvirk, Pichon & Teyssier 2008; Kereš et al. 2009). GB13 therefore implicitly assume that all of the satellite galaxies have been acquired at about that time. This is in conflict with the findings by Gao et al. (2004), according to which about 90 per cent of the present-day subhaloes in cosmological simulations have been accreted after z = 1. The survival mass fraction for haloes accreted at z = 2 is just 2 per cent. Tidal effects would have strongly affected these satellites, in particular if they were moving on radial orbits. This would make the satellites pass through the centre of M31, where they would be tidally disrupted and where they would suffer orbital decay due to dynamical friction (Angus, Diaferio & Kroupa 2011). The scenario presented by GB13 therefore requires the ad hoc assumption that satellites accreted via cold streams are much more stable than usual dark matter subhaloes.

  • The analysis of GB13 implies that the M31 satellites share a few radial orbits, because they were accreted radially from only three to eight different directions. However, they do not test whether this is consistent with the observational data. We have done so using two different methods which are given below.

    • The scenario of GB13 implies that the M31 satellites all originate from no more than eight different directions, as this is the maximum number of streams in their analysis. Because the satellites are supposed to move on radial orbits, their positions as seen from M31 today should unambiguously identify their orbits (which in this scenario are the lines connecting the centre of M31 with the satellites). We can also deduce a satellite's orbital phase from its observed line-of-sight velocity relative to M31's line-of-sight velocity. If a satellite on a radial orbit is closer to the MW than M31 and is approaching us faster than M31, it has to be receding from M31. Thus, its infall direction has to be on the opposite side of M31, so in an all-sky-plot of M31 the mirrored position of the satellite identifies the infall direction. Deducing the original infall direction in this way for all galaxies within 500 kpc of M31, we arrive at the all-sky-plot of putative infall directions shown in Fig. 7. If the assumptions of GB13 were correct, the infall directions of the satellite galaxies should cluster very close to a small number of common directions. This is not the case, demonstrating that the satellite galaxies do not share a few radial orbits, and therefore cannot have been accreted along a few radial streams.

    • As a second test we can make use of the simplicity of the orbital geometry implied by GB13. Assuming a radial orbit, a satellite's measured line-of-sight velocity (in the M31 rest frame) is just the projection of its radial velocity towards (or away from) M31 on to the line of sight. Therefore, we can determine the radial velocity of the satellite from its position relative to M31 and its line-of-sight velocity. As the uncertain satellite position relative to M31 has significant influence on the resulting radial velocity, we use a Monte Carlo method to vary the positions of the satellite and of M31 within their uncertainties. The resulting absolute radial velocities |vrad| of the M31 satellites are plotted against their distances dM31 from M31 in Fig. 8. Some of the satellite galaxies, both within the GPoA and outside, would have very high absolute velocities of more than 300 km s−1 under that assumption. They are likely exceeding the escape velocities from point masses of 1.0 and 2.0 × 1012 M, marked by the two red dotted lines.7 Such satellites would therefore not have remained bound to M31 if they were accreted on radial orbits at z ≥ 2 as considered by GB13, in particular as the galaxy's mass was probably considerably smaller at that time. The fact that the satellites are still close to M31 today therefore requires them to move on non-radial orbits, such that they are in conflict with the basic assumption of GB13.

  • Assuming that the MW and M31 satellites originate from cosmological or cold stream accretion, there is no reason to consider that the MW has acquired its satellite galaxies in a completely different way than M31. Of the 11 satellite galaxies for which proper motion measurements are available, all but one have a dominating tangential velocity component (see e.g. fig. 3 in Pawlowski & Kroupa 2013), and they all have radial velocities well below 200 km s−1 (e.g. Fouquet et al. 2012; Barber et al. 2014). The hypothesis of radial orbits is therefore not supported by the observed velocities of the MW satellites.

  • Cold streams are not expected to be perfectly radial. For example, Wang et al. (2014) assume an impact parameter of 30 kpc. The angular momentum of streams brought into a satellite system is non-negligible. Therefore, the stream material and possible satellites accreted along it will share a common orbital plane, but not be on a radial orbit. Furthermore, we have to consider the observed co-orbiting sense of the GPoA, i.e. that 13 out of the 15 satellites have the same orbital sense (Ibata et al. 2013). This has been ignored in the analysis of GB13. The preference of the GPoA satellites to co-orbit is an integral part of the peculiarity of that structure and further emphasizes the importance of incorporating the orbital direction of accreted satellite galaxies in the analysis.

  • At least one cold stream will feed the formation of the central galactic disc. Gas clumps brought in with the stream merge with the central galaxy and the stream will not be available to bring in satellite galaxies (Dekel et al. 2009).

  • Cosmological zoom-in simulations of the cosmic web (e.g. Ocvirk et al. 2008; Ceverino, Dekel & Bournaud 2010) from which the cold accretion scenario of galaxy mass assembly has been implied (Dekel et al. 2009) show that clumpy substructures develop within the gas filaments. Depending on the spatial and mass resolutions of the simulations, the clump masses can be understood as representing high-mass satellites. GB13 follow this interpretation. The separation of these density structures in simulations varies between 25 and 100 kpc. If one such cold filamentary stream should have caused the accretion of the satellite system of M31, the accumulation of 50 substructures (less than twice the number of the observed satellites to allow for the accretion of clumps by the parent galaxy) separated by 50 kpc with an infall velocity of 200 km s−1 would take about 12.5 Gyr. Since the cold stream infall would be preferably directed towards the LG's centre of mass, the MW and M31 would move along their orbits and by this change the impact parameter and direction of the cold streams, invalidating the assumptions of GB13. The accretion of all satellites along one or a small number of streams is therefore unrealistic.

Figure 7.

All-sky plot as seen from the centre of M31. The symbols indicate the infall directions of M31 satellites determined from their positions relative to M31 and their line-of-sight velocities, assuming them to be on perfectly radial orbits. The circles indicate satellites which belong to the GPoA as identified by Ibata et al. (2013), stars belong to the other satellites in the Ibata et al. (2013) sample and crosses indicate the remaining known galaxies within 500 kpc from M31. The lines indicate the position uncertainties of the infall directions corresponding to the 1σ uncertainties in the heliocentric distances. The red dotted sections of lines indicate that the satellite would not be infalling to but receding from M31 in that direction, because its distance uncertainty makes it change from the side closer to the MW than M31 to the far side or vice versa. While several of the satellites belonging to the GPoA have similar infall directions, close to the position of the MW indicated by the green diamond, most of the satellite galaxy infall directions are widely separated. This is in contrast with the scenario assumed in GB13, in which they all share a few radial orbits. The assumption of all satellite being accreted radially along a few cold streams is not valid if one intends to explain the GPoA with cold stream accretion.

Figure 8.

Median absolute velocities of the M31 satellite galaxies assuming that they are on radial orbits with respect to M31, plotted against their median distances from M31. The error bars indicate the 1σ quantiles of the respective quantities derived from varying the satellite and M31 positions within their observational uncertainties. Symbols are the same as in Fig. 7. Several galaxies have velocities that are likely exceeding the escape velocities of point masses of 1.0 and 2.0 × 1012 M, marked as the dotted red lines. This makes it very unlikely that all satellites have been accreted along a few common cold streams at z ≥ 2, as the unbound ones would have left the vicinity of M31 in the meantime. For comparison it shall be mentioned that all but one of the MW satellites for which proper motions have been estimated have dominating tangential velocity components and that they all have radial velocities well below 200 km s−1 (e.g. Fouquet et al. 2012; Pawlowski & Kroupa 2013; Barber et al. 2014).

Many of these limitations and inconsistencies with the observed situation are due to the oversimplified assumption that the satellite galaxies are confined to a few perfectly radial orbits, which GB13 have modelled as straight lines. A number of these limitations have already been pointed out by GB13, but they did not account for them in their analysis. In the following, we modify their analysis by assuming that each stream defines an orbital plane instead of a radial orbital line. This also allows us to make use of the observational information regarding the co-orbiting satellites in the GPoA.

4.2 New method

We adopt the hypotheses of GB13 that (i) all satellites are accreted via cold streams, (ii) 3–8 streams are feeding the host halo, (iii) the streams are randomly oriented and (iv) each stream is loaded with an equal amount of satellites. Our analysis still makes the (unreasonable) assumption that the origin of the streams does not change, which would spread out the planes. This assumption therefore overestimates the likelihood to find a sufficiently thin satellite plane. As discussed before it is in conflict with cosmological simulations and dynamical friction to assume that all present-day satellites have been accreted before z = 2, when cold streams might have been an important accretion mechanism. We will therefore discuss in Section 4.3 how the situation changes when assuming that half of the satellites have been accreted from directions unrelated to the streams.

Instead of modelling the cold streams as lines crossing the centre of the host galaxy, which assumes that all satellite galaxies are accreted along and remain on perfectly straight radial orbits, it is better to describe the streams as orbital planes. A stream flowing not exactly towards a host galaxy but having a (small) impact parameter will result in a common orbital plane for the stream material. The exact impact parameter is not important in this description. It only has to be assumed that the streams will stay within their original orbital plane. Precession effects in non-spherically symmetric potentials would invalidate this assumption and result in a thickening of the satellite planes (Klarmann et al. in preparation), such that this assumption again results in an overestimation of the probability to find sufficiently thin satellite planes. In our analysis, each stream is described by its pole, which is the direction of the angular momentum of the stream material. This stream pole is a normal vector to the stream's orbital plane, but also contains information about the orbital sense of accreted streams.8

For each realization, we generate M randomly oriented stream poles by drawing them from an isotropic distribution, following the recipe of GB13. We then run through all possible combinations of n = [2, 3, …, M] streams and measure the largest angle between all the stream poles in each of the |$\frac{M!}{(M-n)!\,n!}$| possible combinations for each of the M − 1 different subsample sizes n. For each combination of streams we determine the maximum angle between the stream poles. This is the opening angle of that stream combination. We then record the smallest opening angle, ψ, among all stream combinations of a given subsample size n. This is equivalent to the opening angle of the best-aligned group of streams. As we use the stream poles and not axial normal vectors, this assures that the streams are co-orbiting if ψ < 90°. For each number of streams M = [3, 4, …, 7, 8] we generate 105 random realizations.

For the GPoA, GB13 use an opening angle of 3|$_{.}^{\circ}$|6. We double this value in our analysis. Thus, all subsamples with ψGPoA ≤ 7|$_{.}^{\circ}$|2 are considered to result in planes as narrow as the observed GPoA. From the 27 satellite galaxies in the sample of Ibata et al. (2013), 15 are aligned in the GPoA but only 13 of these co-orbit. As we, in contrast to the analysis of GB13, consider the orbital sense of the accreted satellites, we have to find nGPoA = 13/27 of the M planes to be aligned to ψ ≤ ψGPoA to reproduce this observed situation (which assumes that all streams bring in the same number of satellites and ignores the two aligned but counter-orbiting galaxies in the GPoA). In determining the appropriate nGPoA for a given number of M streams we furthermore have to round up to the next integer, i.e. nGPoA(M) = ⌈13/27 × M⌉. Rounding down would not result in a sufficient number of satellites within the GPoA, while rounding up allows some of the satellites to have been scattered out of the plane.

4.3 Results

Fig. 9 compiles the results of our analysis. Each panel shows the cumulative distribution function of the best-aligned n stream orbit planes for a different number of streams M plotted against the opening angle ψ of the best-aligned combination of planes for a total of 105 realizations.

Figure 9.

Cumulative probabilities to find at least n out of M randomly oriented stream poles within an angle ψ or less of each other. Each panel presents the results of a different value for M, the total number of streams feeding the host. In each panel, the cumulative probabilities of the different n are colour coded, the number nGPoA(M) = ⌈13/27 × M⌉ corresponding to the observed one is marked as a thicker line. The dashed line indicates the adopted opening angle of the observed GPoA around M31, ψGPoA = 7|$_{.}^{\circ}$|2. In contrast to the findings of GB13, who assumed that all M31 satellites move along lines passing through M31's centre, our results demonstrate that it is very unlikely to find ≥13/27 of the satellite galaxies within a thin, co-orbiting plane consistent with the observed one if the streams are described as orbital planes.

A situation in which a sufficient number nGPoA of streams aligns in a co-orbiting plane with an opening angle smaller than or equal to that of the GPoA is very unlikely (see the probabilities listed in Table 6). The largest probability PGPoA = 2.2 per cent that nGPoA (M = 4) = 2 streams align to at least 7|$_{.}^{\circ}$|2 is found for M = 4 streams.

Table 6.

Likelihoods to find a sufficient number of aligned stream orbital planes.

M nGPoA PGPoA [per cent] P(n = nGPoA − 1) [per cent]
3 2 1.157 100.0
4 2 2.227 100.0
5 3 0.006 3.874
6 3 0.015 5.830
7 4 0.000 0.024
8 5 0.000 0.035
M nGPoA PGPoA [per cent] P(n = nGPoA − 1) [per cent]
3 2 1.157 100.0
4 2 2.227 100.0
5 3 0.006 3.874
6 3 0.015 5.830
7 4 0.000 0.024
8 5 0.000 0.035

Probabilities PGPoA to find nGPoA out of M randomly oriented stream planes to be aligned in a co-orbiting orientation. P(n = nGPoA − 1) gives the probability to find nGPoA − 1 aligned planes, which however is not sufficient to explain the observed number of co-orbiting satellites in the GPoA.

Table 6.

Likelihoods to find a sufficient number of aligned stream orbital planes.

M nGPoA PGPoA [per cent] P(n = nGPoA − 1) [per cent]
3 2 1.157 100.0
4 2 2.227 100.0
5 3 0.006 3.874
6 3 0.015 5.830
7 4 0.000 0.024
8 5 0.000 0.035
M nGPoA PGPoA [per cent] P(n = nGPoA − 1) [per cent]
3 2 1.157 100.0
4 2 2.227 100.0
5 3 0.006 3.874
6 3 0.015 5.830
7 4 0.000 0.024
8 5 0.000 0.035

Probabilities PGPoA to find nGPoA out of M randomly oriented stream planes to be aligned in a co-orbiting orientation. P(n = nGPoA − 1) gives the probability to find nGPoA − 1 aligned planes, which however is not sufficient to explain the observed number of co-orbiting satellites in the GPoA.

Overall, an alignment of a sufficient number of streams is essentially excluded for M ≥ 5 as the probabilities PGPoA remain below 0.015 per cent in these cases. However, according to simulations it is most likely that galaxies are fed by a lower number of streams. Nevertheless, a situation comparable to the observed GPoA occurs in only 1 out of 50 galaxies, if assuming that all satellites are accreted exclusively in streams (see Section 4.1 for why this is unrealistic).

While the likelihoods to have fewer streams than necessary to align are higher (right-hand column in Table 6), this is insufficient to reproduce the number of co-orbiting satellites in the GPoA, as we have to require that at least 13 out of 27 satellites are in a common plane. Satellites can be scattered out of the plane later (as was also argued by GB13), but it is less likely to add objects to the co-orbiting plane. Thus it is incorrect to include cases where the fraction of aligned planes is smaller than nGPoA(M) = ⌈13/27 × M⌉.

As discussed in Section 4.1 it is expected that additional satellite galaxies have been accreted after the cold stream era. The addition of such satellites that fall in unrelated to the cold streams requires more objects to be in the common plane initially if the GPoA is to be explained by the cold stream scenario. If we assume for simplicity that 50 per cent of all satellites (0.5 × 27 = 13.5) have been accreted individually9 and that these follow (to first order) an isotropic distribution. Then only half of the M31 satellites would have been accreted in streams. This in turn requires that essentially all stream planes had to align initially. The reason is that, if drawn from isotropy, only sin (3|$_{.}^{\circ}$|6) = 6.3 per cent of isotropically distributed satellites will end up within the GPoA, i.e. about 1 out of 13.5. As we also consider the corotating trend of the satellites, this number has to be divided by 2 once more, because there is only a 50 per cent chance that a satellite with a random direction of motion follows the common velocity trend. Thus, we can assume that on average 0.5 satellites out of the 13 co-orbiting in the GPoA are from the random distribution, while the remaining 12.5 have to originate from the alignment of cold stream orbital planes. As we assume that 50 per cent of the satellites (= 13.5) are accreted via cold streams, at least 12.5 out of these 13.5 had to align initially. Thus, essentially all cold stream planes had to be aligned. As can be seen in Fig. 9, the probability for this to happen is negligible for all assumed numbers of streams.

4.4 Summary

Given these results, we disagree with the conclusions of GB13. We have shown that with the simple but well-motivated change from radial to planar satellite orbits we are already unable to reproduce their results. A sufficient alignment of streams occurs rarely (in at most ≈2 per cent of the systems with four streams, less in all other configurations). This finding stands in contrast to the results by GB13 who report probabilities as high as 50 to 100 per cent.

5 CONCLUSION

We have searched for co-orbiting planes of satellite galaxies similar to the GPoA around M31 and the VPOS around the MW within satellite systems based on the Millennium II cosmological simulation. Our analysis was motivated by and is similar to the works by BB14 and W13. However, we have identified several limitations in their sample selection and methods, which we have addressed in our analysis.

Using the same simulation data, we find different results than BB14 concerning the frequency of the satellite plane around M31. A GPoA-like plane, i.e. a satellite structure which is sufficiently narrow, radially extended and contains at least 13 satellite galaxies with the same orbital sense, is found with a probability of only PGPoA = 0.04 to 0.17 per cent, depending on the adopted criteria and assumed distance to Andromeda XXVII. This has to be contrasted to the 2 per cent found by BB14, who only required two of the three properties to be similar to the observed situation. Our analysis of their data therefore disputes their conclusion, according to which ‘[…] Millennium II haloes with 13 or more corotating satellites orbiting in a thin plane are not rare […]’ and ‘[…] the existence of a [GPoA] is not in conflict with the standard cosmological framework’. The GPoA is extremely unlikely to be found in a dark matter based ΛCDM universe within the limitations of current semi-analytic galaxy formation models. Its existence around M31 therefore poses a serious challenge to that and similar cosmological models of galaxy formation. These conclusions are supported and further strengthened by the analysis of the Millennium II simulation by Ibata et al. (2014), who apply the exact PAndAS survey area to the simulated satellite systems and require the satellite planes to be seen in the same edge-on orientation as the observed GPoA. Using this very different method, they find that satellite planes sharing the same parameters as the GPoA are extremely unlikely among simulated satellite systems (0.04 per cent).

We have also searched for satellite planes similar to the VPOS around the MW in the same data set, incorporating the effect of randomly oriented obscuring MW discs. This analysis is comparable to the one conducted by W13. We are unable to reproduce their results on the frequency of satellite systems similarly flattened as observed, finding a 2 to 3 times lower frequency using only one criterion. To be comparable with the observed VPOS, a satellite system furthermore has to be not only as flattened as the 11 classical MW satellites, but also sufficiently radially extended, and coherently orbiting within the best-fitting plane. Only a fraction PVPOS ≤ 0.06 per cent of the simulated satellite systems fulfil these requirements simultaneously. In contrast to the conclusion by W13, our results evidence that the observed phase-space correlated satellite galaxy system of the MW alone already poses a serious challenge to the ΛCDM model, if it is well represented by the analysed Millennium II simulation in combination with semi-analytic modelling.

Another attempt to reconcile the observed GPoA with cosmological expectations has been made by GB13, who argue that the accretion of satellite galaxies along cold gas streams naturally results in planar satellite distributions. Their model of satellites moving on perfectly radial orbits along cold streams is an oversimplification, is in conflict with the observed M31 satellite galaxies, and is unable to produce the predominantly co-orbiting GPoA. We have adjusted their analysis by requiring all satellites accreted along a cold stream to remain in the same orbital plane, instead of on perfectly radial orbits, i.e. straight lines. With this modification we are unable to reproduce their results because it is much less likely that a sufficient number of randomly oriented planes align to the degree required to be consistent with the narrow, co-orbiting GPoA. Assuming further that only half of the satellite galaxies have been accreted along streams results in a negligible likelihood to obtain a structure as thin as the GPoA.

Both satellite galaxy systems for which we have sufficient observational data, that of the MW and of M31, are therefore in conflict with ΛCDM expectations derived from dark-matter-only simulations combined with a semi-analytic galaxy formation model. The reason for this can in principle either indicate a failure of the semi-analytic models on dwarf galaxy scale, missing physical processes in the simulation, or of the underlying cosmology (see also Kroupa 2012).

The GPoA and the VPOS are both similarly unlikely to be found in the same satellite systems based on the Millennium II simulation, such that both structures can be considered to be of similar significance. As we have based the GPoA and the VPOS analyses on the same set of host galaxies, we are now also able to determine the ΛCDM likelihood that both systems contain such structures.

Formally, taking the most favourable value for the GPoA of PGPoA = 0.0017 the probability of having a pair of two galaxies displaying one GPoA and one VPOS, in any order, is P = 2 × PGPoA × PVPOS = 2.04 × 10−6.

Nevertheless, one could argue that any individual configuration of satellite systems of two galaxies, if looked at with too much scrutiny, will have almost zero probability. To avoid this and be more conservative, we can consider that, since the two satellite planes have different properties (most importantly the GPoA consists of only a subset of the M31 satellites whereas almost all MW satellites are situated in the VPOS), the likelihood to have either a GPoA or a VPOS around any one random host galaxy is the sum of the individual probabilities:10
\begin{equation*} P_{\mathrm{sum}} \le P_{\mathrm{VPOS}} + P_{\mathrm{GPoA}} = 2.3 \times 10^{-3}. \nonumber \end{equation*}
This situation happens precisely twice in the two major galaxies of the LG, which has a probability
\begin{equation*} P_{\mathrm{LG}} = P_{\mathrm{sum}}^2 \le 5.29 \times 10^{-6}. \nonumber \end{equation*}
The chance to find the VPOS- and/or GPoA-like satellite planes in the LG is therefore extremely small if all satellites are primordial dwarf galaxies in a ΛCDM universe.

Our work challenges recent claims that the observed satellite galaxy structures are consistent with the ΛCDM model. Co-orbiting planes of satellite galaxies such as the GPoA around M31 and the VPOS around the MW are still in conflict with satellite galaxy distributions derived from ΛCDM simulations. In fact, our results indicate that due to improved knowledge of the phase-space correlation of the systems they are even more in conflict with those than before.

We would like to thank Henning Bahl and Holger Baumgardt for providing the Millennium II data used in their study and for providing additional information on their analysis. We furthermore thank them for discussing their initial results on the MW and M31 satellite systems with MSP and PK in Summer 2013. We also thank Tobias Goerdt for discussing his work with GH. MSP, PK, FL and HJ acknowledge financial support from DAAD/Go8 travel grant 56265912. HJ also acknowledges the financial support from the Australian Research Council through the Discovery Project grant DP120100475. DM was supported by the National Science Foundation under grant no. AST 1211602 and by the National Aeronautics and Space Administration under grant no. NNX13AG92G. JD acknowledges funding through FONDECYT grant no. 3140146. DAF acknowledges the financial support from the Australian Research Council through the Discovery Project grant DP130100388.

1

The feature which we call the VPOS is called the ‘Magellanic plane’ by Kunkel & Demers (1976) and Lynden-Bell (1976); the ‘great pancake’ by Libeskind et al. (2005), Li & Helmi (2008) and Wang, Frenk & Cooper (2013); and the disc of satellites by Kroupa, Theis & Boily (2005) and Metz, Kroupa & Libeskind (2008).

2

The same structure has been called the vast thin disc of satellites by Hammer et al. (2013) and the vast thin plane of dwarf galaxies by Bahl & Baumgardt (2014).

3

Note that the notion of ‘rare is already here slightly subjective, as a 2 per cent occurrence could certainly already be labelled as ‘not very common’.

4

However, as we stick with the exact same simulation data set, we do not address the potential effect of selecting galaxies residing in clusters rather than completely isolated groups.

5

BB14 do not mention the adopted value, but looking at their fig. 2 it appears to be approximately 190 kpc.

6

Even though that work in its present form is only a non-refereed arXiv preprint, it is being cited (including by published refereed papers) as indicative that cold stream accretion might explain the GPoA (e.g. Shaya & Tully 2013, BB14; Hoffmann et al. 2014).

7

Treating M31 as a point mass is a conservative approximation as well, as for a given total mass, the mass distribution with the highest escape velocity is a point mass.

8

In fact we would have to require at least one stream to feed the growth of M31's galactic disc, but refrain from incorporating this possibility for the sake of simplicity. Later merger events might also change the host galaxy's spin direction.

9

Note that this is a very low fraction given that, according to Gao et al. (2004), 90 per cent of the subhaloes that survived until today have been accreted after z = 1, whereas simulations indicate that cold stream are dominant only before a redshift of z = 2.

10

Some of the VPOS-like satellite systems might simultaneously be GPoA-like, thus the sum is only an upper limit.

REFERENCES

Angus
G. W.
Diaferio
A.
Kroupa
P.
MNRAS
2011
, vol. 
416
 pg. 
1401
 
Bahl
H.
Baumgardt
H.
MNRAS
2014
, vol. 
438
 pg. 
2916
  
(BB14)
Barber
C.
Starkenburg
E.
Navarro
J. F.
McConnachie
A. W.
Fattahi
A.
MNRAS
2014
, vol. 
437
 pg. 
959
 
Barnes
J. E.
Hernquist
L.
Nature
1992
, vol. 
360
 pg. 
715
 
Bellazzini
M.
Oosterloo
T.
Fraternali
F.
Beccari
G.
A&A
2013
, vol. 
559
 pg. 
L11
 
Bett
P.
Eke
V.
Frenk
C. S.
Jenkins
A.
Okamoto
T.
MNRAS
2010
, vol. 
404
 pg. 
1137
 
Bournaud
F.
, et al. 
Science
2007
, vol. 
316
 pg. 
1166
 
Bournaud
F.
Duc
P.-A.
Emsellem
E.
MNRAS
2008
, vol. 
389
 pg. 
L8
 
Boylan-Kolchin
M.
Springel
V.
White
S. D. M.
Jenkins
A.
Lemson
G.
MNRAS
2009
, vol. 
398
 pg. 
1150
 
Brainerd
T. G.
ApJ
2005
, vol. 
628
 pg. 
L101
 
Ceverino
D.
Dekel
A.
Bournaud
F.
MNRAS
2010
, vol. 
404
 pg. 
2151
 
Conn
A. R.
, et al. 
ApJ
2012
, vol. 
758
 pg. 
11
 
Conn
A. R.
, et al. 
ApJ
2013
, vol. 
766
 pg. 
120
 
Cooper
A. P.
, et al. 
MNRAS
2010
, vol. 
406
 pg. 
744
 
D'Onghia
E.
Lake
G.
ApJ
2008
, vol. 
686
 pg. 
L61
 
Dabringhausen
J.
Kroupa
P.
MNRAS
2013
, vol. 
429
 pg. 
1858
 
Deason
A. J.
, et al. 
MNRAS
2011
, vol. 
415
 pg. 
2607
 
Dekel
A.
Birnboim
Y.
MNRAS
2006
, vol. 
368
 pg. 
2
 
Dekel
A.
, et al. 
Nature
2009
, vol. 
457
 pg. 
451
 
Duc
P.-A.
, et al. 
MNRAS
2011
, vol. 
417
 pg. 
863
 
Duc
P.-A.
Paudel
S.
McDermid
R. M.
Cuillandre
J.-C.
Serra
P.
Bournaud
F.
Cappellari
M.
Emsellem
E.
MNRAS
2014
, vol. 
440
 pg. 
1458
 
Forbes
D. A.
Kroupa
P.
Metz
M.
Spitler
L.
Mercury
2009
, vol. 
38
 pg. 
24
 
Fouquet
S.
Hammer
F.
Yang
Y.
Puech
M.
Flores
H.
MNRAS
2012
, vol. 
427
 pg. 
1769
 
Galianni
P.
Patat
F.
Higdon
J. L.
Mieske
S.
Kroupa
P.
A&A
2010
, vol. 
521
 pg. 
A20
 
Gao
L.
White
S. D. M.
Jenkins
A.
Stoehr
F.
Springel
V.
MNRAS
2004
, vol. 
355
 pg. 
819
 
Goerdt
T.
Burkert
A.
2013
 
arXiv:1307.2102v1 (GB13)
Guo
Q.
, et al. 
MNRAS
2011
, vol. 
413
 pg. 
101
 
Guo
Q.
White
S.
Angulo
R. E.
Henriques
B.
Lemson
G.
Boylan-Kolchin
M.
Thomas
P.
Short
C.
MNRAS
2013
, vol. 
428
 pg. 
1351
 
Hammer
F.
Yang
Y. B.
Wang
J. L.
Puech
M.
Flores
H.
Fouquet
S.
ApJ
2010
, vol. 
725
 pg. 
542
 
Hammer
F.
Yang
Y.
Fouquet
S.
Pawlowski
M. S.
Kroupa
P.
Puech
M.
Flores
H.
Wang
J.
MNRAS
2013
, vol. 
431
 pg. 
3543
 
Hoffmann
K.
, et al. 
2014
 
Ibata
R. A.
, et al. 
Nature
2013
, vol. 
493
 pg. 
62
 
Ibata
R. A.
Ibata
N. G.
Lewis
G. F.
Martin
N. F.
Conn
A.
Elahi
P.
Arias
V.
Fernando
N.
ApJ
2014
, vol. 
784
 pg. 
L6
 
Jarosik
N.
, et al. 
ApJS
2011
, vol. 
192
 pg. 
14
 
Karachentsev
I. D.
Bautzmann
D.
Neyer
F.
Polzl
R.
Riepe
P.
Zilch
T.
Mattern
B.
2014
 
Keller
S. C.
Mackey
D.
Da Costa
G. S.
ApJ
2012
, vol. 
744
 pg. 
57
 
Kereš
D.
Katz
N.
Fardal
M.
Davé
R.
Weinberg
D. H.
MNRAS
2009
, vol. 
395
 pg. 
160
 
Koch
A.
Grebel
E. K.
AJ
2006
, vol. 
131
 pg. 
1405
 
Kroupa
P.
PASA
2012
, vol. 
29
 pg. 
395
 
Kroupa
P.
Theis
C.
Boily
C. M.
A&A
2005
, vol. 
431
 pg. 
517
 
Kroupa
P.
, et al. 
A&A
2010
, vol. 
523
 pg. 
A32
 
Kunkel
W. E.
Demers
S.
Dickens
R. J.
Perry
J. E.
Smith
F. G.
King
I. R.
R. Greenwich Obser. Bull., Vol. 182
1976
The Magellanic Plane
pg. 
241
 
Lemson
G.
Virgo
Consortium
2006
 
astro-ph/0608019
Libeskind
N. I.
Frenk
C. S.
Cole
S.
Helly
J. C.
Jenkins
A.
Navarro
J. F.
Power
C.
MNRAS
2005
, vol. 
363
 pg. 
146
 
Libeskind
N. I.
Frenk
C. S.
Cole
S.
Jenkins
A.
Helly
J. C.
MNRAS
2009
, vol. 
399
 pg. 
550
 
Libeskind
N. I.
Knebe
A.
Hoffman
Y.
Gottlöber
S.
Yepes
G.
Steinmetz
M.
MNRAS
2011
, vol. 
411
 pg. 
1525
 
Li
Y.-S.
Helmi
A.
MNRAS
2008
, vol. 
385
 pg. 
1365
 
Lovell
M. R.
Eke
V. R.
Frenk
C. S.
Jenkins
A.
MNRAS
2011
, vol. 
413
 pg. 
3013
 
Lynden-Bell
D.
MNRAS
1976
, vol. 
174
 pg. 
695
 
McConnachie
A. W.
AJ
2012
, vol. 
144
 pg. 
4
 
McConnachie
A. W.
Irwin
M. J.
MNRAS
2006
, vol. 
365
 pg. 
902
 
Metz
M.
Kroupa
P.
MNRAS
2007
, vol. 
376
 pg. 
387
 
Metz
M.
Kroupa
P.
Jerjen
H.
MNRAS
2007
, vol. 
374
 pg. 
1125
 
Metz
M.
Kroupa
P.
Libeskind
N. I.
ApJ
2008
, vol. 
680
 pg. 
287
 
Metz
M.
Kroupa
P.
Theis
C.
Hensler
G.
Jerjen
H.
ApJ
2009
, vol. 
697
 pg. 
269
 
Nelson
D.
Vogelsberger
M.
Genel
S.
Sijacki
D.
Kereš
D.
Springel
V.
Hernquist
L.
MNRAS
2013
, vol. 
429
 pg. 
3353
 
Ocvirk
P.
Pichon
C.
Teyssier
R.
MNRAS
2008
, vol. 
390
 pg. 
1326
 
Paudel
S.
, et al. 
ApJ
2013
, vol. 
767
 pg. 
133
 
Pawlowski
M. S.
Kroupa
P.
MNRAS
2013
, vol. 
435
 pg. 
2116
 
Pawlowski
M. S.
McGaugh
S. S.
MNRAS
2014
, vol. 
440
 pg. 
908
 
Pawlowski
M. S.
Kroupa
P.
de Boer
K. S.
A&A
2011
, vol. 
532
 pg. 
A118
 
Pawlowski
M. S.
Pflamm-Altenburg
J.
Kroupa
P.
MNRAS
2012a
, vol. 
423
 pg. 
1109
 
Pawlowski
M. S.
Kroupa
P.
Angus
G.
de Boer
K. S.
Famaey
B.
Hensler
G.
MNRAS
2012b
, vol. 
424
 pg. 
80
 
Pawlowski
M. S.
Kroupa
P.
Jerjen
H.
MNRAS
2013
, vol. 
435
 pg. 
1928
 
Ploeckinger
S.
Hensler
G.
Recchi
S.
Mitchell
N.
Kroupa
P.
MNRAS
2014
, vol. 
437
 pg. 
3980
 
Recchi
S.
Theis
C.
Kroupa
P.
Hensler
G.
A&A
2007
, vol. 
470
 pg. 
L5
 
Richardson
J. C.
, et al. 
ApJ
2011
, vol. 
732
 pg. 
76
 
Sales
L.
Lambas
D. G.
MNRAS
2004
, vol. 
348
 pg. 
1236
 
Sales
L.
Lambas
D. G.
MNRAS
2009
, vol. 
395
 pg. 
1184
 
Sawa
T.
Fujimoto
M.
PASJ
2005
, vol. 
57
 pg. 
429
 
Shaya
E. J.
Tully
R. B.
MNRAS
2013
, vol. 
436
 pg. 
2096
 
Springel
V.
, et al. 
MNRAS
2008
, vol. 
391
 pg. 
1685
 
Wang
J.
Frenk
C. S.
Cooper
A. P.
MNRAS
2013
, vol. 
429
 pg. 
1502
  
(W13)
Wang
L.
Zhu
W.
Feng
L.-L.
Macciò
A. V.
Chang
J.
Kang
X.
MNRAS
2014
, vol. 
439
 pg. 
L85
 
Wetzstein
M.
Naab
T.
Burkert
A.
MNRAS
2007
, vol. 
375
 pg. 
805
 
Yang
Y.
Hammer
F.
ApJ
2010
, vol. 
725
 pg. 
L24
 
Zhao
H.
Famaey
B.
Lüghausen
F.
Kroupa
P.
A&A
2013
, vol. 
557
 pg. 
L3