The dawn – dusk asymmetry of ion density in the dayside magnetosheath and its annual variability measured by THEMIS

The local and global plasma properties in the magnetosheath play a fundamental role in regulating solar wind–magnetosphere coupling processes. However, the magnetosheath is a complex region to characterise as it has been shown theoretically, observationally and through simulations that plasma properties are inhomogeneous, non-isotropic and asymmetric about the Sun-Earth line. To complicate matters, dawn–dusk asymmetries are sensitive to various changes in the upstream conditions on an array of timescales. The present paper focuses exclusively on dawn–dusk asymmetries, in particularly that of ion density. We present a statistical study using THEMIS data of the dawn–dusk asymmetry of ion density in the dayside magnetosheath and its long-term variations between 2009 and 2015. Our data suggest that, in general, the dawn-side densities are higher, and the asymmetry grows from noon towards the terminator. This trend was only observed close to the magnetopause and not in the central magnetosheath. In addition, between 2009 and 2015, the largest asymmetry occurred around 2009 decreasing thereafter. We also concluded that no single parameter such as the Alfvén Mach number, plasma velocity, or the interplanetary magnetic field strength could exclusively account for the observed asymmetry. Interestingly, the dependence on Alfvén Mach number differed between data sets from different time periods. The asymmetry obtained in the THEMIS data set is consistent with previous studies, but the solar cycle dependence was opposite to an analysis based on IMP-8 data. We discuss the physical mechanisms for this asymmetry and its temporal variation. We also put the current results into context with the existing literature in order to relate THEMIS era measurements to those made during earlier solar cycles.


Introduction
The complex and nonlinear interaction between the solar wind and the terrestrial magnetosphere presents a challenging and rich geophysical problem. The goal is to fully understand the redistribution of the solar wind kinetic energy as it undergoes many transitions before its impact on our immediate geospace environment and the ionosphere. Firstly, the terrestrial magnetic field presents an abrupt obstacle to the incoming solar wind flow causing it to rapidly slow. This deceleration forms a fast mode shock wave that stands upstream of the planet and effectively re-distributes the solar wind kinetic energy into other degrees of freedom (Sagdeev, 1966;Sagdeev and Galeev, 1969;Papadopoulos, 1985). Immediately downstream of the bow shock, the magnetosheath transition layer houses the "shocked" solar wind plasma which has been significantly altered by the bow shock front. This processing means that in general, the magnetosheath plasma is hotter, denser, slower, more turbulent and magnetically reconfigured. Therefore it is no surprise that the parametrisation of the upstream-downstream transition has been the focus of investigations over the past several decades using a combination of theoretical models, numerical simulations and in situ observations. However, despite the many advances made during this time, and considering the relatively close proximity between the upstream and downstream plasma, many open questions remain. One unresolved problem of the magnetosheath is related to the presence of numerous dawn-dusk asymmetries. To summarise these asymmetries: magnetic field strength, plasma velocity, temperature anisotropy and mirror mode occurrence favour the dusk-flank and magnetic field turbulence, ion total temperature and ion number density are stronger on the dawn-side (Paularena et al., 2001;Němeček et al., 2003;Longmore et al., 2005;Walsh et al., 2012;Dimmock and Nykyri, 2013;Dimmock et al., 2014Dimmock et al., , 2015a. Remaining issues are the following:

Published by Copernicus
(1) what are their underlying physical mechanisms, and (2) what role do they play in driving magnetospheric plasma properties such as the formation of the cold dense plasma sheet. The latter point relates directly to plasma transport processes other than those which are reconnection-driven. These processes are more relevant during prolonged periods of northward interplanetary magnetic field (IMF) when subsolar reconnection (Dungey, 1961) has moved to higher latitudes (Dungey, 1963). To complicate matters, the means of transporting mass, momentum and energy across the magnetopause involves a synthesis of magnetohydrodynamic (MHD), kinetic, and numerous wave-particle interaction processes. In spite of that, several candidates capable of facilitating transport across the magnetopause have been identified: magnetic reconnection (Dungey, 1961), impulsive penetration (Karlsson et al., 2012;Archer and Horbury, 2013;Plaschke et al., 2013), Kinetic Alfvén waves (Johnson and Cheng, 1997;Chaston et al., 2008) the Kelvin Helmholtz Instability (KHI) (Miura and Pritchett, 1982;Miura, 1992;Chen et al., 1997;Otto and Fairfield, 2000;Nykyri and Otto, 2001;Nykyri et al., 2006;Hasegawa et al., 2009) and flux transfer events (Russell and Elphic, 1978). It is noteworthy that each of these processes do not strictly operate independently, and during certain conditions, plasma transport may be a result of a combination of these processes operating simultaneously. Of particular interest with regard to the present paper is that the efficiency of plasma transport processes like the KHI depend on the global and local magnetosheath conditions (Nykyri, 2013). Therefore, it is likely that magnetosheath dawn-dusk asymmetries play a substantial role in driving magnetospheric conditions mainly by (1) asymmetries inherent to the natural seed population, and (2) by regulating plasma entry processes such as KHI.
Recently, Nykyri (2013) presented compelling evidence that the dawn-dusk asymmetry of the magnetic field strength has a direct impact on the growth of the KHI at the night-side magnetopause. Using global and local MHD simulations, the author showed that when the magnetic field strength tangential to the Kelvin-Helmholtz k-vector was weaker, the Alfvén speed was reduced, which increased the velocity shear ratio (V shear /V A × k). As a result, during a Parker-spiral IMF, the dawn-side magnetopause should be more favourable to the growth of the KHI since the (steady state) tangential magnetic field strength was weaker on the quasi-parallel magnetosheath flank. This result also took into account the faster magnetosheath flow (and enhanced velocity shear) on the dusk-flank (Longmore et al., 2005;Walsh et al., 2012;Dimmock and Nykyri, 2013) during these conditions. Interestingly, Taylor et al. (2012) reported Cluster observations showing that the occurrence of the KHI was larger on the dusk-flank when observations were limited to the dayside. The results of Nykyri (2013) are relevant close to the terminator and on the nightside. Although a direct comparison between the studies is difficult, an intriguing spatial dependence is suggested. Recently, Walsh et al. (2015) presented evidence that dense plasma on the magnetospheric side originating from a plasmaspheric plume can cause the magnetopause to be more unstable to the KHI during cases when the velocity shear is low. The authors argued that under the plume conditions, the instability might develop under lower velocity shears allowing the KHI to form more sunward than usual. What is clear from these studies is that the development of the KHI is sensitive to the local variations of plasma properties adjacent to the magnetopause. For that reason, it is important to quantify the dawn-dusk asymmetries, their dependence on upstream conditions, as well as their short-or long-term temporal variability. Paularena et al. (2001) used IMP eight measurements to derive the dawn-dusk asymmetry of ion density over the range −20 < X GSE < −15 R E . A strong dawn-favoured asymmetry was found but appeared to obey a complicated temporal dependency. The authors investigated two periods from 1978 to 1980 and between 1994 to 1997 which corresponded to solar maximum and minimum, respectively. The largest asymmetry was identified during solar maximum whereas little to no asymmetry was recorded during solar minimum. Interestingly, the authors found no clear evidence that the IMF orientation, solar wind speed, or Alfvén Mach number variations could be causing the differences between the two data sets. Němeček et al. (2003) investigated the dawn-dusk asymmetry of magnetosheath ion fluxes derived from Interball-1 at more Sunward locations between −15 and +5 R E . The authors also found a strong dawnfavoured asymmetry, and concluded that it could be related, but not strictly driven by, the IMF orientation. Longmore et al. (2005) using Cluster data also ruled out the IMF as the cause for a dawn-favoured density asymmetry. Using a combination of Time History of Events and Macroscale Interactions during Substorms (THEMIS) dayside measurements just outside the magnetopause, and Block-Adaptive-Tree-Solarwind-Roe-Upwind-Scheme (BATS-R-US) MHD simulations, Walsh et al. (2012) reported a dawn-favoured asymmetry of approximately 20 % close to the terminator. Interestingly, since Walsh et al. (2012) selected data close to the solar minimum during 2008-2010, the 20 % asymmetry reported appears contradictory to the solar cycle dependence found during the previous solar minimum by Paularena et al. (2001). However, the discrepancy in the results could also suggest a difference in the conditions on the dayside and nightside. Dimmock and Nykyri (2013)  not detect any strong dawn-dusk asymmetry in their analysis of the THEMIS data, which could again indicate a temporal dependence as their study spanned the solar minimum from October 2007 through the rising phase until December 2012. Another explanation could be that the asymmetry was too delicate to be identified from their processing of the data. To conclude, it is clear that open questions remain regarding the dawn-dusk asymmetry of the ion density, and that the effects of the temporal variations of magnetosheath ion densities on the KHI and the associated plasma transport are not fully understood.
The present study aims to expand upon previous work by utilising the full THEMIS database (2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015) to quantitatively study the dawn-dusk asymmetry of magnetosheath ion density. We utilised our existing data analysis tool developed over the past 3 years (see Dimmock and Nykyri, 2013;Dimmock et al., 2014Dimmock et al., , 2015a to compile annual data sets of magnetosheath ion density measurements between 2008 through 2014. We also investigate the solar wind parameter dependence for each year separately to shed further light on the physical driver. We discuss the findings in terms of the physical mechanisms, and compare the results with earlier findings.
The manuscript is structured as follows. Section 2 describes the data sources and the methods adopted in processing these data. The final statistical data sets are presented in Sect. 3 in which we describe the trends observed in each one. We discuss the physical aspects of these results in Sect. 4 and put them into context with the existing literature. Finally, we draw our conclusions and provide a brief summary in Sect. 5. Some additional information is provided in Appendices A and B which are supplementary to the results in Sect. 3.

Data sets and instrumentation
Our magnetosheath data set were compiled entirely from measurements conducted between October 2007 and December 2014 using the suite of instrumentation onboard each of the THEMIS spacecraft (Angelopoulos, 2008). For this study, ion densities were obtained through the onboard moments provided by the ElectroStatic Analyser (ESA) instrument . We used the level-2 (L2) data files which provide measurements at ∼ 3 s intervals corresponding to the spacecraft spin period. For the filtering of our data set, and for the input parameters to the magnetosheath boundary models, we extracted the solar wind parameters from NASA/GSFC's OMNI data set through OM-NIWeb (http://omniweb.gsfc.nasa.gov). We used the highresolution (1 min cadence) observations that have been propagated (King and Papitashvili, 2005) from the measuring spacecraft locations to the bow shock nose (Farris and Russell, 1994). To account for density changes in the magnetosheath caused by density differences in the upstream solar wind, all density values were normalised by their corresponding upstream value so that n i = n ms /n sw . All results shown are based on the normalised database unless otherwise stated.

Compilation of magnetosheath density measurements
For each THEMIS probe location, we apply a transformation from the Geocentric Solar Ecliptic (GSE) frame to the Magnetosheath InterPlanetary Medium (MIPM) reference frame (Bieber and Stone, 1979;Verigin et al., 2006;Dimmock and Nykyri, 2013). The purpose of this transformation is to automatically arrange the data in a manner which permits a direct comparison of observations during vastly different solar wind conditions. The MIPM frame arranges data based on the upstream IMF vector such that the bow shock geometry of the magnetosheath dawn and dusk flanks remain quasiparallel and quasi-perpendicular, respectively. The MIPM magnetosheath closely resembles the magnetosheath configuration during a typical Parker-spiral IMF. For example, if a magnetosheath observation during an ortho-Parker-spiral IMF was made on the dawn-flank, this would in the MIPM frame be located on the dusk-side. Therefore, to reflect this arrangement, from this point we adopt the notation of and ⊥ to describe data collected on the quasi-parallel (dawn) and quasi-perpendicular (dusk) flanks, respectively. The IMF orientation (for each magnetosheath point) is assessed based on a 20 min average of the OMNI time series data. The 20 min window length was selected to account first, the unknown convection time to each probe location, second, the inaccuracies associated with the solar wind measurements, and third, to minimise the impact from erroneously timed transient solar wind events. To account for the motion of the magnetosheath boundaries, the fractional distance F MIPM across the model magnetosheath is measured for each point, where |R| corresponds to the radial spacecraft location, r mp and r bs are the radial distances to the model magnetopause and bow shock positions (along R), respectively. The values of F MIPM range from 0 at the magnetopause, to 1 at the bow shock, such that the intermediate values indicate the relative position between these boundaries. The instantaneous location of the magnetosheath boundaries were calculated using the Shue et al. (1998) magnetopause model, and the Verigin et al. (2001Verigin et al. ( , 2003Verigin et al. ( , 2006 bow shock model based on input variables provided by the 20 min averaged OMNI data. Each data-point is also corrected for planetary aberration by orienting the x axis (in each MIPM rotation) along the solar wind flow velocity vector with the planetary orbital motion subtracted from the y component. The magnetosheath plasma parameters were obtained by evaluating a 3-min centre-weighted average at each location. In this case, the mean value of the ion number density is computed.

A. P. Dimmock et al.: Magnetosheath ion density asymmetry
Within a given 3 min window, if values further than the standard deviation from the mean σ/ n i > 1 were present, these windows were excluded. As a precaution, we also removed time instances when the magnetosheath flow exceeded that of the solar wind. This criteria is set for the sole purpose of eliminating solar wind outliers close to the bow shock arising from erroneous bow shock location identification. In the current context, this criterion removed less than 0.1 % of our data set; and thus, has no bearing on the results. Regarding magnetospheric points, we remove points which lie very far from their typical distribution and are consistent with magnetospheric values. Since we do not have simultaneous measurements in the magnetosphere, we employ no criterion to eliminate magnetospheric points which could be close to those observed in the magnetosheath. However, we take many steps to validate our results and to ensure that they are not a manifestation of magnetospheric outliers; these are discussed later in the manuscript. For each datapoint, the MIPM location, simultaneous solar wind conditions, and the magnetosheath plasma parameters were stored to create a database of magnetosheath measurements spanning October 2007 through December 2014. Temporal or parameter-based subsets could then be extracted from this complete database.
In this case, we isolated yearly subsets to monitor the annual variations of the ion density dawn-dusk asymmetry.
The initial 1-year window (1 January-31 December 2008) was then incremented every 6 months (1 July 2008-1 July 2009) until the final window spanned 1 January-31 December 2014. As a result, each window (except the first and last) contained a 50 % overlap between the previous and next windows. The purpose of the overlap was to provide improved temporal resolution over the complete THEMIS interval and to minimise edge effects from window averaging. Since the spatial coverage of the THEMIS probes in the magnetosheath were not consistent over the entire period (see Fig. A1 in Appendix A), some physical restrictions were enforced to ensure that the yearly subsets contained comparable spatial distributions of the measurements. We restricted our data collection to the dayside magnetosheath, discarding data around the subsolar region (0 • ) and at the terminator (90 • ). The extent of this exclusion is a 15 • sector at both locations. As a result, our statistical data cover angles between 15 • and 75 • with respect to local noon and the terminator. In addition, we only apply data collected close to the magnetopause defined by the following range of fractional distances, 0 ≤ F MIPM ≤ 0.2. Figure 1 shows an illustration of the data collection area highlighted by the orange shaded regions, showing also the angular sectors excluded from the analysis. This region was chosen, as it has continuous coverage throughout the observing period as indicated by the yearly distribution of THEMIS measurements ( Fig. A1 in Appendix A). Furthermore, we targeted the region close to the magnetopause as our previous results (Dimmock and Nykyri, 2013) indicated that while no clear asymmetry was present in the central dayside magnetosheath (0.33 ≤ F MIPM ≤ 0.66), the data suggested an asymmetry could be present closer to the magnetopause. Walsh et al. (2012) also utilised THEMIS measurements within close proximity to the magnetopause and reported a dawn-favoured asymmetry. In addition, the plasma properties in close proximity to the magnetopause are the ones regulating plasma transport processes.
To measure the dawn-dusk asymmetry of each annual subset, the mean and median of ion densities on the dawn (n ) and dusk (n ⊥ ) flanks were computed within the confines of the limits defined above. Because of the large and uneven variability of the upstream conditions, we cannot further divide the data sets to smaller angular bins, as the bins then would not have similar solar wind statistics. For the dawn and dusk flank data sets, we check the solar wind statistics by plotting the probability density functions of the corresponding solar wind parameters for each dawn and dusk flank (see Figs. B1-B3) of each yearly subset to ensure they are comparable between the data sets. The averaged densities are then used as indicators of the density on each flank. The asymmetry (A) is then computed as follows Therefore, A in Eq.
(2) represents the percentage of which the average dawn density is larger than the dusk. The nature of the asymmetry can be inferred from the polarity of A such that The asymmetry is computed for each annual data set to provide an estimate of the asymmetry as a function of the year A(t).

Solar wind measurements
Studies aimed at characterising the interaction between the solar wind and planetary bodies are complicated by the dynamic nature of the processes across a large range of temporal and spatial scales. For example, the polarity of the IMF B z component has a profound impact through magnetic reconnection recorded as responses in geomagnetic indices (Osmane et al., 2015) over time periods much less than 1 hour. On the other hand, the analysis of long-term trends present in spacecraft observations suggest that much longer (∼ 11 year) variations originating from changes in the solar activity also produce a measurable effect (e.g., Paularena et al., 2001).  It should be mentioned that σ was computed from the 20 min averages and therefore excludes information about the highfrequency variations. The following results characterise the driving of the ion density dawn dusk-asymmetry based on the solar wind input parameters during the THEMIS observation period shown as the shaded blue area. It is clear even from this 7-year interval that the solar wind parameters also vary over slow (yearly) timescales; which may be reflected in the magnetosheath conditions. presents the dawn-dusk asymmetry computed as described in Eqs. (2) and (3). Since the data coverage is sufficient, we computed the asymmetry as a function of the angle from local noon A(θ ) and no angular bins were excluded. Please note that angles of 0 and 90 • correspond to local noon and the terminator, respectively. The grey bars show the mean estimates whereas the orange bars are calculated from the median. The angular resolution of each bar is 15 • which is then incremented by a 50 % overlap to form the next bin. As mentioned previously, the data in panels (b and c) include data only within 0.2F MIPM of the magnetopause. The blue error-bars were calculated from the maximum possible variation of the asymmetry resulting from the standard error of the mean. Panel (c) shows the mean and median crosssectional cuts of the data used to estimate the asymmetry presented in panel (b). The error-bars here also originate from the standard error within each bin. In general, the mean and median averages yielded the same result, and the asymmetries provided are comfortably within the calculated bounds of error. There is a dawn-favoured asymmetry which ranges from around 1 % at local noon to approximately 19 % at the terminator. In general, the asymmetry grows with increasing angular distance from the subsolar region. The average asymmetry measured between local noon and the terminator is approximately 8 and 10 % for mean and median averages, respectively. Based on data collected in the central magnetosheath (not shown), the asymmetry alternated between the dawn and dusk flanks, and so these data proved inconclusive. Our data are therefore suggestive that the dawn-dusk asymmetry is dependent upon the spatial location in which it is measured throughout the magnetosheath. This spatial dependence is both a function of the radial distance across the magnetosheath and the distance from the subsolar region, i.e. A(F MIPM , θ ). Figure 4 shows the dawn-dusk asymmetry of ion density between 2008 though 2014, A(t 2008→2014 ). The grey and orange bars in Fig. 4a indicate the annual asymmetry based on mean and median averages, respectively. For comparative purposes, the solid and dashed blue lines show the asymmetry calculated over the entire database. Each error-bar is calculated in the same manner as that in Fig. 3b. A time series plot of solar wind parameters for M A , |V|, |B| and n i are shown underneath in panel (b) as indicated by the grey, red, green and blue solid lines, respectively. These solar wind parameters are 30-day averages of the entire OMNI data set and not limited to the subset corresponding to the magnetosheath datapoints. We also checked the corresponding upstream conditions for this subset, and they were very similar to the full data set. What is striking in Fig. 4a is the large asymmetry recorded from 2008 throughout mid-2010. These years encompass the 2009 solar minimum, evidenced by the values of the solar wind properties plotted in Fig. 4b below. The corresponding solar wind parameters during 2008-2010 are the following: M A ∼ 11, |V| ∼ 400 km s −1 , |B| ∼ 3 nT and n i ∼ 5 cm −3 . During the period 2011-2014, the asymmetry reduced to less than 5 %. An increase during the period 2012-2013 does not reach the values recorded by THEMIS during the solar minimum interval. Because the data presented here are collected through an automated process in close proximity to the magnetopause, it is important to ensure that the asymmetry is not artificially created by outliers introduced by the inaccuracy of the magnetopause model. Since the magnetopause is not an infinitesimally thin boundary and its position is a complex function of the solar wind and magnetospheric dynamics, outliers of magnetospheric origin close to the magnetopause can be expected. For this reason, we plot the probability density function of each yearly data set in Fig. 5 which clearly demonstrates that the core of each dawn-flank distribution is shifted to higher density values than its dusk-side counterpart. In addition, this also demonstrates that the asymmetry values (based on lower order moments) do not result from comparing ill-matched distribution functions between each flank.

Dawn-dusk asymmetry Alfvén Mach number dependence
The data presented in Fig. 6 describe the dawn-dusk asymmetry as a function of the upstream Alfvén Mach number, A(M A ). Since M A changes from ∼ 11 during 2009 to ∼ 8.5 in 2013 (see Figs. 2 and 4), then it is logical to examine its role in regulating the asymmetry. The grey bars in Fig. 6 correspond to measurements during 2008-2009 whereas the orange bars represent data during 2010-2014. These intervals differ slightly from the ones before, however these provided adequate data for each interval while still splitting them between solar minimum and the rising phase. The error-bars are calculated in the same manner as explained previously.
For reference, the average dawn-dusk asymmetry for both intervals are indicated by the solid and dashed green horizontal lines. The vertical orange lines mark the average M A for both intervals. Both panels show the same data but the asymmetries were calculated using mean (Fig. 2a) and median ( Fig. 2b) averages, respectively. Each bar corresponds to an average M A using a window length 3. This window is then incremented by 1.5 from 5 to 15. As before, overlapping windows allow better resolution over the parameter range. The total (dawn + dusk) amount of data present in the 2008-2009 and 2010-2014 data sets for each M A window were [11029, 17955, 23023, 21802, 16631, 11 467, 6705] and [24 501, 31 991, 30 122, 22 786, 16 308, 9929, 6383] points, respectively for M A = [6. 5, 8.0, 9.5, 11.0, 12.5, 14.0, 15.5]. The data were relatively evenly distributed between the dawn and dusk flanks. What is immediately obvious from Fig. 6 is that around solar minimum, the asymmetry appears to be inversely correlated with M A . What is puzzling is that for the rising solar maximum period (2011)(2012)(2013)(2014)(2015) there appears to be no clear Mach number dependence, and the asymmetry is relatively unchanged over the entire M A regime. These dependencies are true of both the mean and median processed data. We varied the widths of the M A bins and obtained the same result. A similar observation was also reported by Paularena et al. (2001) in which the M A dependence differed between their temporal data sets; we will return to this point in the discussion.

Discussion
The aim of this investigation was to accurately quantify the dawn-dusk asymmetry of ion density in the dayside magnetosheath and determine its time dependence using THEMIS measurements. To compile our data sets, we extracted magnetosheath observations from the complete THEMIS catalogue using our magnetospheric data analysis tool (see Dimmock and Nykyri (2013); Dimmock et al. (2014Dimmock et al. ( , 2015a). We quantitatively showed the asymmetry during typical solar wind conditions (Fig. 3), which was shown to be more pronounced and clearer closer to the magnetopause compared to the central magnetosheath. To study the annual dependence of the asymmetry, we segmented the complete data set into 12-month subsets for which we computed the asymmetry. To compensate for the varying data coverage of the THEMIS probes over the time period, we confined our statistical data collection to a fixed region such that comparable data coverage was reached between the data sets. All data were collected on the dayside magnetosheath for angles of 15 • ≤ θ ≤ 75 • , with respect to local noon (0 • ) and the terminator (90 • ). In addition, we only included data collected close to the magnetopause within 20 % of the magnetosheath thickness (fractional distance range of 0 ≤ F MIPM ≤ 0.2, with respect to the magnetopause (0) and the bow shock (1); see Fig. 1 and Fig. A1 in Appendix A). The dawn-dusk asymmetry was then evaluated (see Eqs. 2 and 3) for each subset to determine the temporal variability (Fig. 4). As the earlier results on the asymmetry have been partially conflicting, it is worth mentioning that the results are statistically significant and that the distribution functions on each flank have statistically similar properties (Fig. 5). We also confirmed that the asymmetries do not arise from different solar wind conditions during measurements made in the dawn and dusk flanks (see Figs. B2 and B3) or from inaccuracies in the magnetopause determination. If a mismatch of outliers were present in the data sets, the mean and median would be skewed by longer tails generated by the very low densities even if the distribution cores would overlap; this was shown not to be the case. In addition, we also ensured that the results were not driven by the size of the data sets. To investigate this, we randomly sampled each subset so they were of equal size. These data yielded remarkably similar results, indicating that the orbital changes over the mission duration did not affect our results. As a final check, we systematically Ann. Geophys., 34, 511-528 varied the limits of F MIPM (from 0.1 to 0.3) and θ , effectively changing the data selection region shown by the orange region in Fig. 1 by moving away from the magnetopause. We found that even though we eliminated data very close to the magnetopause (0 < F MIPM < 0.1) our conclusions were unchanged; suggesting magnetospheric data were not responsible for driving the asymmetry. It is also noteworthy that our estimates for 2008-2010 are remarkably similar to that of Walsh et al. (2012) using the same data set but different methodology; suggestive of a robust result. We showed that the magnetosheath ion density is higher on the dawn flank near the magnetopause (Fig. 3b). We also examined the central magnetosheath, but in this region no clear trends were identified. These results are qualitatively and quantitatively similar to those by Walsh et al. (2012) who also examined the asymmetry using THEMIS data close to the magnetopause. Both studies also suggest that the asymmetry grows as a function of distance from local noon, in concert (but not necessarily driven) with the increase in asymmetrical magnetosheath thickness. Our results are in agreement with Paularena et al. (2001) and Longmore et al. (2005), who also reported higher densities on the dawnside magnetosheath. The repeatability of the asymmetry using several spacecraft (THEMIS, Cluster, and IMP 8) and different processing methodologies strongly suggests that the asymmetry is of physical origin as opposed to a manifestation from statistical or technical means. It is also worth reiterating that the MIPM magnetosheath is similar to the GSE magnetosheath during a Parker-spiral IMF; therefore these results should be statistically relevant as this is the most common IMF orientation.
Over the THEMIS era, the maximum asymmetry was observed during and around the solar minimum year 2009. From 2011 onwards, the asymmetry reduced and became more ambiguous. To further examine the solar cycle dependence, we investigated the role played by the solar wind Alfvén Mach number (see Figs. 2 and 6), but these data did not suggest a clear link between the two. Around solar minimum, the maximum asymmetry (22.5 %) occurred during lower Alfvén Mach number (M A = 6-7) while a much lower asymmetry (2.5 %) was observed for higher Alfvén Mach numbers (M A = 13-14). After 2010, there were no clear trends, and the asymmetry varied less than 5 percentage points over the entire regime of Alfvén Mach numbers. While the Alfvén Mach number likely plays a role in regulating the asymmetry to some extent, there also seems to be other factors contributing to its long-term variability. The Alfvén Mach number dependence can either arise from the modification of the downstream plasma properties as indicated by the Rankine-Hugoniot jump conditions (Kulsrud, 2005), or by the asymmetrical changes in the magnetosheath thickness (Wu, 1992). The asymmetrical magnetosheath thickness alone has been shown to produce dawndusk asymmetries in MHD simulations using the BATS-R-US code (Walsh et al., 2012;Dimmock and Nykyri, 2013).
One further thing to note is that around solar minimum, the standard deviation of the Alfvén Mach number was enhanced (see Fig. 2 iii-3), but the possible impact of this effect on the asymmetry is unclear.
It is well established that solar wind-magnetosphere coupling undergoes seasonal (semiannual) variations (Russell and McPherron, 1973) and therefore, by means of feedback to the magnetosheath, could have an impact on magnetosheath conditions close to the magnetopause. This is particularly important since the THEMIS probes sampling of the dawn and dusk flanks (in the GSE & GSM frames) are seasonally dependent. If these feedback processes were strong, then even though solar wind data were comparable over 12months intervals, an asymmetry could manifest in the data set. However, we do not believe these are responsible for the trends we observed. We computed the dawn-dusk asymmetry of the dawn-dusk sampling, and by comparing this with the ion density asymmetry variations, this did not offer any conclusive or clear relation. In addition, this is a good time to re-iterate the arrangement of data in the MIPM frame. Even during periods when THEMIS occupied only a single flank, data would still be placed on the appropriate flank based on the IMF orientation. Therefore, in general, both flanks are sampled over the entire year (albeit not equally) which should reduce such seasonal effects. Finally, since the magnetosphere is best described in a GSM frame, and we adopt the MIPM frame in the magnetosheath, matching or quantifying such feedback processes is not straightforward, and in the scope of the current investigation, would be speculative at best. However, it is possible they may play a role, and further work is needed to identify and quantify their statistical impact on magnetosheath properties.
Temporal variation of the ion density dawn-dusk asymmetry was previously reported by Paularena et al. (2001) using magnetosheath observations from IMP 8, together with solar wind observations from ISEE 1, ISEE 3, and WIND. It should be noted that the data analysed by Paularena et al. (2001) were collected between −20 ≤ X GSE ≤ −15 R E while our study focused entirely on the dayside. In addition, the orbit of IMP 8 covered much higher latitudes in which Z GSE varied between ±20 R E whereas THEMIS is closer to the equatorial plane. In spite of these differences, both studies present strong evidence of a dawn-favoured asymmetry, which exhibit significant variations as a function of the solar cycle phase. However, the nature of the temporal dependence is not consistent between the two investigations; Paularena et al. (2001) covered the solar maximum during August 1978-February 1980 and the solar minimum during November 1994 to October 1997. In contrast to our results indicating the strongest asymmetry during solar minimum, Paularena et al. (2001) obtained a maximum asymmetry during the solar maximum time period. In agreement with our results, Paularena et al. (2001) report a clearer dependence on Mach number during solar minimum, and no discernible relationship to upstream parameters during solar maximum. The authors also conclude that the Alfvén Mach number was not responsible for driving the asymmetry differences between the solar cycle maximum and minimum. While the final reason for the contradicting results of the asymmetry remains unclear, a partial explanation could be that the plasma conditions during the two minima and maxima periods were markedly different. For example, dynamic pressure, density, velocity and magnetic field strength are visibly lower during 2009 as opposed to 1996 (see Figs. 2 and 2), and the values during the 2011-2014 maximum period are at the same level as during the minimum of 1996. Thus, although the ion density dawn-dusk asymmetry does vary with the 11-year solar cycle, the behaviour of the plasma properties during each solar cycle have to be taken into account. With the increased sophistication and resolution of numerical models and simulations (Alfthan et al., 2014;Karimabadi et al., 2014;Palmroth et al., 2015), in conjunction with the long-term monitoring of solar wind conditions, we should be well equipped to resolve the drivers of this asymmetry.
To summarise: this paper shows strong evidence of a dawn-favoured asymmetry of ion density throughout the dayside magnetosheath and that its magnitude varies as a function of time. The physical processes leading to magnetosheath dawn-dusk asymmetries (Paularena et al., 2001;Longmore et al., 2005;Walsh et al., 2012;Dimmock and Nykyri, 2013;Dimmock et al., 2014Dimmock et al., , 2015a, are fundamental in solar wind-magnetosphere coupling as they can directly impact plasma transport processes at the magnetopause. Our results imply that the role they play is time dependent and this needs to be taken into consideration in the analysis of long duration data sets of both magnetosheath and magnetospheric measurements. Furthermore, the present results are of importance to numerical simulations (Alfthan et al., 2014;Karimabadi et al., 2014;Palmroth et al., 2015) attempting to accurately reproduce the magnetosheath behaviour in response to upstream conditions. Open questions remain on the detailed solar wind driver of the ion density asymmetry, and further work is required to use the state-ofthe-art numerical models to resolve the causal relationships.

Conclusions
We have studied the dawn-dusk asymmetry of ion density in the dayside magnetosheath and quantified its annual variability using THEMIS data. Our conclusions can be summarised as follows: 1. Ion densities were higher on the dawn-side magnetosheath close to the magnetopause (average ∼ 8 %) but this trend was unclear in the central magnetosheath.
2. The asymmetry increased with angular distance from the subsolar point and peaked at the terminator (∼ 20 %).
3. Our results are consistent with previous studies reporting a dawn-favoured asymmetry (Longmore et al., 2005;Walsh et al., 2012), but the solar cycle dependence in our data is different from that reported by Paularena et al. (2001).
4. In our data set, the peak asymmetry (∼ 20 %) was observed during 2009 coinciding with the solar minimum.
5. The asymmetry notably decreased and became inconclusive in the years 2011-2014 which took place after the 2009 minima.
6. We examined the role of the Alfvén Mach number and found a larger asymmetry during low values at solar minimum, but for the years 2011-2014 we identified no clear trend.
We conclude that additional work is required to explain the solar cycle relationship measured by different data sets. In addition, these results should be of consequence to the regulation of viscous plasma transport processes (such as KHI), since they are sensitive to the local plasma properties in close vicinity to the magnetopause. Finally, magnetosheath dawndusk asymmetries and their temporal dependencies should be kept in mind for future investigations which include such processes, and in particular when they utilise observations over a long period of time.

Data availability
The OMNI data were accessed by NASA/GSFC's Space Physics Data Facility's OMNIWeb. These data can be accessed at http://omniweb.gsfc.nasa.gov. THEMIS data were accessed via http://themis.ssl.berkeley.edu/index.shtml. All data are available free of charge. The colour in each statistical map corresponds to the number of points within each 0.5 R E square bin. Please note, these data are plotted in a base 10 logarithmic scale. As discussed in the main body of the manuscript, the 60 • angular sector in which the data were used are indicated in each panel. Figure A1 shows the yearly coverage of THEMIS magnetosheath measurements in the MIPM frame. The colour scale represents the number of points per 0.5 × 0.5 bin.  Figure B3. Yearly probability density functions of solar magnetic field strength (|B|) which correspond to magnetosheath data points on the dawn and dusk flanks. The "2" in each yearly label indicates that the data set spans the middle of that year to the middle of the next. e.g. 2012 2 = July 2012-July 2013.