The Milky Way Stellar Halo Is Twisted and Doubly Broken: Insights from DESI DR2 Milky Way Survey Observation Songting Li1,2,3aa, Wenting Wang1,3aa, Sergey E. Koposov4,5aa, João A. S. Amarante1,3aa, Alis J. Deason6aa, Nathan R. Sandford7,8aa, Ting S. Li7,8aa, Gustavo E. Medina7,8aa, Jiaxin Han1,3aa, Monica Valluri9aa, Oleg Y. Gnedin9aa, Namitha Kizhuprakkat10, Andrew P. Cooper10aa, Leandro Beraldo e Silva11aa, Carlos Frenk6aa, Raymond G. Carlberg7,8aa, Mika Lambert12aa, Tian Qiu1,3aa, Jessica Nicole Aguilar13aa, Steven Ahlen14aa, Davide Bianchi15,16aa, David Brooks17aa, Todd Claybaugh13, Axel de la Macorra18aa, Peter Doel17aa, Jaime E. Forero-Romero19,20aa, Enrique Gaztañaga21,22aa, Satya Gontcho A Gontcho13,23aa, Gaston Gutierrez24aa, Dick Joyce25, Robert Kehoe26aa, Anthony Kremin13aa, Claire Lamman27aa, Martin Landriau13aa, Laurent Le Guillou28aa, Ramon Miquel29,30aa, Will Percival31,32,33aa, Francisco Prada34aa, Ignasi Pérez-Ràfols35aa, Graziano Rossi36, Eusebio Sanchez37aa, David Schlegel13aa, Ray Sharples6,38aa, Joseph Harry Silber13aa, David Sprayberry25aa, Gregory Tarlé39aa, Benjamin Alan Weaver25, and Hu Zou40aa 1 Department of Astronomy, School of Physics and Astronomy, and Key Laboratory for Particle Astrophysics and Cosmology (MOE)/Shanghai Key Laboratory for Particle Physics and Cosmology, Shanghai Jiao Tong University, Shanghai 200240, People’s Republic of China; songtingli@sjtu.edu.cn, wenting.wang@sjtu.edu.cn 2 Tsung-Dao Lee Institute, Shanghai Jiao Tong University, Shanghai 201210, People’s Republic of China 3 State Key Laboratory of Dark Matter Physics, School of Physics and Astronomy, Shanghai Jiao Tong University, Shanghai 200240, People’s Republic of China 4 Institute for Astronomy, University of Edinburgh, Royal Observatory, Blackford Hill, Edinburgh EH9 3HJ, UK; Sergey.Koposov@ed.ac.uk 5 Institute of Astronomy, University of Cambridge, Madingley Road, Cambridge CB3 0HA, UK 6 Institute for Computational Cosmology, Department of Physics, Durham University, South Road, Durham DH1 3LE, UK 7 Department of Astronomy & Astrophysics, University of Toronto, Toronto, ON M5S 3H4, Canada 8 David A. Dunlap Department of Astronomy & Astrophysics, University of Toronto, 50 St George Street, Toronto, ON M5S 3H4, Canada 9 Department of Astronomy and Astrophysics, University of Michigan, Ann Arbor, MI, USA 10 Institute of Astronomy and Department of Physics, National Tsing Hua University, Hsinchu 30013, Taiwan 11 Observatório Nacional, Rio de Janeiro - RJ, 20921-400, Brazil 12 Department of Astronomy & Astrophysics, University of California, Santa Cruz, 1156 High Street, Santa Cruz, CA 95064, USA 13 Lawrence Berkeley National Laboratory, 1 Cyclotron Road, Berkeley, CA 94720, USA 14 Department of Physics, Boston University, 590 Commonwealth Avenue, Boston, MA 02215, USA 15 Dipartimento di Fisica “Aldo Pontremoli,” Università degli Studi di Milano, Via Celoria 16, I-20133 Milano, Italy 16 INAF-Osservatorio Astronomico di Brera, Via Brera 28, 20122 Milano, Italy 17 Department of Physics & Astronomy, University College London, Gower Street, London, WC1E 6BT, UK 18 Instituto de Física, Universidad Nacional Autónoma de México, Circuito de la Investigación Científica, Ciudad Universitaria, Cd. de México C. P. 04510, Mexico 19 Departamento de Física, Universidad de los Andes, Cra. 1 No. 18A-10, Edificio Ip, CP 111711, Bogotá, Colombia 20 Observatorio Astronómico, Universidad de los Andes, Cra. 1 No. 18A-10, Edificio H, CP 111711 Bogotá, Colombia 21 Institut d’Estudis Espacials de Catalunya (IEEC), c/ Esteve Terradas 1, Edifici RDIT, Campus PMT-UPC, 08860 Castelldefels, Spain 22 Institute of Cosmology and Gravitation, University of Portsmouth, Dennis Sciama Building, Portsmouth, PO1 3FX, UK 23 University of Virginia, Department of Astronomy, Charlottesville, VA 22904, USA 24 Fermi National Accelerator Laboratory, PO Box 500, Batavia, IL 60510, USA 25 NSF NOIRLab, 950 N. Cherry Avenue, Tucson, AZ 85719, USA 26 Department of Physics, Southern Methodist University, 3215 Daniel Avenue, Dallas, TX 75275, USA 27 The Ohio State University, Columbus, OH 43210, USA 28 Sorbonne Université, CNRS/IN2P3, Laboratoire de Physique Nucléaire et de Hautes Energies (LPNHE), FR-75005 Paris, France 29 Institut de Física d’Altes Energies (IFAE), The Barcelona Institute of Science and Technology, Edifici Cn, Campus UAB, 08193, Bellaterra (Barcelona), Spain 30 Institució Catalana de Recerca i Estudis Avançats, Passeig de Lluís Companys, 23, 08010 Barcelona, Spain 31 Department of Physics and Astronomy, University of Waterloo, 200 University Avenue W, Waterloo, ON N2L 3G1, Canada 32 Perimeter Institute for Theoretical Physics, 31 Caroline Street N, Waterloo, ON N2L 2Y5, Canada 33 Waterloo Centre for Astrophysics, University of Waterloo, 200 University Avenue W, Waterloo, ON N2L 3G1, Canada 34 Instituto de Astrofísica de Andalucía (CSIC), Glorieta de la Astronomía, s/n, E-18008 Granada, Spain 35 Departament de Física, EEBE, Universitat Politècnica de Catalunya, c/Eduard Maristany 10, 08930 Barcelona, Spain 36 Department of Physics and Astronomy, Sejong University, 209 Neungdong-ro, Gwangjin-gu, Seoul 05006, Republic of Korea 37 CIEMAT, Avenida Complutense 40, E-28040 Madrid, Spain 38 Centre for Advanced Instrumentation, Department of Physics, Durham University, South Road, Durham DH1 3LE, UK 39 University of Michigan, 500 S. State Street, Ann Arbor, MI 48109, USA 40 National Astronomical Observatories, Chinese Academy of Sciences, A20 Datun Road, Chaoyang District, Beijing, 100101, People’s Republic of China Received 2025 December 1; revised 2026 February 1; accepted 2026 February 2; published 2026 February 26 Abstract Using K giants from the second data release (DR2) of the Dark Energy Spectroscopic Instrument (DESI) Milky Way Survey, we measure the shape, orientation, radial profile, and density anisotropies of the Milky Way (MW) stellar halo over 8 kpc< rGC < 200 kpc. We identify a triaxial stellar halo (axis ratio 10:8:7), 43° tilted from the disk, showing two break radii at ∼16 and ∼76 kpc, likely associated with Gaia-Sausage/Enceladus and the Large Magellanic Cloud (LMC), respectively. The inner stellar halo (<30 kpc) is oblate and aligned with the disk, The Astrophysical Journal, 999:108 (19pp), 2026 March 1 https://doi.org/10.3847/1538-4357/ae41b9 © 2026. The Author(s). Published by the American Astronomical Society. aaaaaaa Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. 1 https://orcid.org/0000-0002-6469-8263 https://orcid.org/0000-0002-5762-7571 https://orcid.org/0000-0003-2644-135X https://orcid.org/0000-0002-7662-5475 https://orcid.org/0000-0001-6146-2645 https://orcid.org/0000-0002-7393-3595 https://orcid.org/0000-0002-9110-6163 https://orcid.org/0000-0003-0105-9576 https://orcid.org/0000-0002-8010-6715 https://orcid.org/0000-0002-6257-2341 https://orcid.org/0000-0001-9852-9954 https://orcid.org/0000-0001-8274-158X https://orcid.org/0000-0002-0740-1507 https://orcid.org/0000-0002-2338-716X https://orcid.org/0000-0002-7667-0081 https://orcid.org/0000-0002-2527-8899 https://orcid.org/0000-0002-4900-2088 https://orcid.org/0000-0003-0822-452X https://orcid.org/0000-0001-6098-7247 https://orcid.org/0000-0001-9712-0006 https://orcid.org/0000-0002-8458-5047 https://orcid.org/0000-0002-1769-1640 https://orcid.org/0000-0002-6397-4457 https://orcid.org/0000-0002-2890-3725 https://orcid.org/0000-0001-9632-0815 https://orcid.org/0000-0003-3142-233X https://orcid.org/0000-0003-0825-0517 https://orcid.org/0000-0002-7101-697X https://orcid.org/0000-0001-6356-7424 https://orcid.org/0000-0002-6731-9329 https://orcid.org/0000-0003-1838-8528 https://orcid.org/0000-0001-7178-8868 https://orcid.org/0000-0002-6610-4836 https://orcid.org/0000-0002-0644-5727 https://orcid.org/0000-0001-7145-8674 https://orcid.org/0000-0001-6979-0125 https://orcid.org/0000-0002-9646-8198 https://orcid.org/0000-0002-5042-5088 https://orcid.org/0000-0003-3449-8583 https://orcid.org/0000-0002-3461-0320 https://orcid.org/0000-0001-7583-6441 https://orcid.org/0000-0003-1704-0781 https://orcid.org/0000-0002-6684-3997 mailto:songtingli@sjtu.edu.cn mailto:wenting.wang@sjtu.edu.cn mailto:Sergey.Koposov@ed.ac.uk https://doi.org/10.3847/1538-4357/ae41b9 https://crossmark.crossref.org/dialog/?doi=10.3847/1538-4357/ae41b9&domain=pdf&date_stamp=2026-02-26 https://creativecommons.org/licenses/by/4.0/ whereas the outer stellar halo becomes prolate and perpendicular to the disk, consistent with the vast polar structure of MW satellites. The twisted halo may arise from the disk−halo angular momentum shift triggered by the infall of a massive satellite. The anisotropic density distribution of the stellar halo is also measured, with successful reidentification of the Hercules-Aquila Cloud North/South (HAC-N/S) overdensity and the Virgo overdensity (VOD). Break radii are found at 15 and 30 kpc for VOD and HAC-N/S, respectively. We identify the LMC transient density wake with a break radius at 60 kpc in the Pisces overdensity region. We also find new observational evidence of the LMC collective density wake, by showing a break radius at ∼100 kpc in the northern Galactic cap with a clear density peak at 90 kpc. In the end, we found that more metal-poor halo stars are more radially extended. Our results provide important clues to the assembly and evolution of the MW stellar halo under the standard cosmic structure formation framework. Unified Astronomy Thesaurus concepts: Milky Way stellar halo (1060); Large Magellanic Cloud (903); Galaxy stellar content (621) 1. Introduction The stellar halo of the Milky Way (MW), although composing only ∼2% of the total stellar mass of the MW (A. J. Deason et al. 2019), encodes a unique record of its accretion history. According to the Λ cold dark matter (ΛCDM) model, the MW stellar halo is formed hierarchically, and the assembly history is embedded in the accretion debris and in the form of stellar overdensities and streams (O. J. Eggen et al. 1962; L. Searle & R. Zinn 1978; S. D. M. White & C. S. Frenk 1991; A. Helmi et al. 1999; K. Freeman & J. Bland-Hawthorn 2002; J. S. Bullock & K. V. Johnston 2005; V. Belokurov 2013; J. Bland-Hawthorn & O. Gerhard 2016; K. V. Johnston 2016; A. Helmi 2020), including, for example, the Hercules-Aquila Cloud North/ South (HAC-N/S) overdensity (V. Belokurov et al. 2007; L. L. Watkins et al. 2009; I. T. Simion et al. 2014), the Virgo overdensity (VOD; A. K. Vivas et al. 2001; H. J. Newberg et al. 2002; M. Jurić et al. 2008; V. Chandra et al. 2023a), and the Pisces overdensity (B. Sesar et al. 2007; J. A. Kollmeier et al. 2009; L. L. Watkins et al. 2009; J. D. Nie et al. 2015; V. Chandra et al. 2023a, 2023b; A. Viswanathan et al. 2024). Contemporary wide-field surveys, from the Gaia space mission (Gaia Collaboration et al. 2016) to large spectroscopic surveys such as SEGUE (B. Yanny et al. 2009), APOGEE (S. R. Majewski et al. 2017), GALAH (S. L. Martell et al. 2017), H3 (C. Conroy et al. 2019), LAMOST (X.-Q. Cui et al. 2012), and DESI (DESI Collaboration et al. 2016a, 2016b, 2022; J. Guy et al. 2023; E. F. Schlafly et al. 2023; T. N. Miller et al. 2024; C. Poppett et al. 2024; DESI Collaboration et al. 2024, 2025; DESI Collaboration et al. 2025), now provide detailed chemical measurements and six-dimensional phase- space information for millions of stars, and thus the density distribution of the MW stellar halo can be directly quantified by star counting. Although assembling a whole-sky stellar sample that probes the Galaxy remains challenging, many investigations have charted the density distribution of the MW stellar halo. Various tracers that have precise distance estimates and are luminous enough can be used to measure the density distribution, including main-sequence turnoff (MSTO; H. L. Morrison et al. 2000; E. F. Bell et al. 2008; B. Sesar et al. 2011; B. Pila-Díez et al. 2015; M. Cavieres et al. 2025), blue horizontal branch (BHB; J. Sommer-Larsen 1987; G. W. Preston et al. 1991; O. Y. Gnedin et al. 2010; A. J. Deason et al. 2011; T. Fukushima et al. 2018; G. F. Thomas et al. 2018; T. Fukushima et al. 2019; E. Starkenburg et al. 2019; F. Yu et al. 2024; J. A. S. Amarante et al. 2024), RR Lyrae (M. R. S. Hawkins 1984; G. W. Preston et al. 1991; C. J. Wetterer & J. T. McGraw 1996; A. K. Vivas & R. Zinn 2006; L. L. Watkins et al. 2009; B. Sesar et al. 2010, 2013; L. Faccioli et al. 2014; N. Hernitschek et al. 2018; G. E. Medina et al. 2018; G. Iorio & V. Belokurov 2019; A. Pieres et al. 2020; K. M. Stringer et al. 2021; G. E. Medina et al. 2024, 2025), and red giant branch (X.-X. Xue et al. 2015; Y. Xu et al. 2018; J. J. Han et al. 2022; C. Yang et al. 2022; J. M. M. Lane et al. 2023; M. López-Corredoira et al. 2024) stars. Most studies model the stellar halo as an oblate spheroid aligned with the Galactic disk, fitting its density with a broken power-law profile (A. J. Deason et al. 2011; C. Yang et al. 2022; J. A. S. Amarante et al. 2024; M. Cavieres et al. 2025). In particular, A. J. Deason et al. (2013) revealed that break radii record apocentric pileups from early major mergers. N. Garavito-Camargo et al. (2019) and R. P. Naidu et al. (2021) further pointed out that major mergers will also induce stellar overdensities in the sky, leading to a highly anisotropic stellar halo. Recently, multiple works (J. J. Han et al. 2022; C. Yang et al. 2022; J. M. M. Lane et al. 2023) have independently recovered a double-break power-law density profile for the MW stellar halo, with the two break radii corresponding to two apocenter passages of Gaia-Sausage/Enceladus (GSE). These discoveries confirmed the dominant role of the GSE merger in shaping the inner ∼30 kpc of the MW stellar halo (V. Belokurov et al. 2018; M. Haywood et al. 2018; A. Helmi et al. 2018; L. Lancaster et al. 2019; R. P. Naidu et al. 2021; W. Wu et al. 2022; J. J. Han et al. 2023; J. Nibauer & A. Bonaca 2025; G. E. Medina et al. 2025). However, the constraints on the exact locations of break radii and the morphology of the MW stellar halo are still subject to large inconsistencies. Based on the LAMOST survey, C. Yang et al. (2022) reported an oblate stellar halo without tilt, and the flattening of the stellar halo is changing with radius. However, J. J. Han et al. (2022; H3) and J. M. M. Lane et al. (2023; APOGEE) have found a prolate stellar halo with tilt. Due to the limited survey depth and sample size, the measurements are mostly within ∼50 kpc in these studies. More recently, and benefiting from the usage of large photo- metric samples to identify BHB fractions, J. A. S. Amarante et al. (2024) studied the anisotropic density distributions of the MW stellar halo by fitting the density distribution at different line-of- sight directions. They also reported that in the outer stellar halo ∼55% of stars reside in ∼20% of sky patches. M. Cavieres et al. (2025) further measured the density profile in these dense sky patches. They detected a clear overdensity in the outer stellar halo, with a break radius at ∼70 kpc. However, both J. A. S. Amarante et al. (2024) and M. Cavieres et al. (2025) use photometric samples, which might suffer from more contaminations and imprecise distance measurements than 2 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. http://astrothesaurus.org/uat/1060 http://astrothesaurus.org/uat/903 http://astrothesaurus.org/uat/621 http://astrothesaurus.org/uat/621 spectroscopic samples (J. J. Han et al. 2022; C. Yang et al. 2022; J. M. M. Lane et al. 2023). The Milky Way Survey (MWS) of DESI will release ∼12 million stellar spectra in its second data release (DR2). It provides one of the largest and deepest spectroscopic samples of the MW stellar halo spanning a wide footprint, and hence it offers an invaluable spectroscopic dataset to investigate the anisotropic density distribution of the MW stellar halo down to unprecedented details. In this paper, we measure the shape, orientation, and density profile of the MW stellar halo out to ∼200 kpc, using the large halo K giant sample from DESI DR2, which is currently an internal version of data release. In particular, we explore how the shape and orientation of the MW stellar halo change with distances, as well as its alignment with the MW disk and the vast polar structures (VPOSs) of MW satellite galaxies, which provides important probes into the assembly and evolution of the MW stellar halo under the standard hierarchical cosmic structure formation scenario. We pay particular attention to the anisotropic density profiles along different directions, look at the imprints from the perturbations of the Large Magellanic Cloud (LMC), and investigate the dependence of the stellar halo density profiles on [Fe/H]. This paper is organized as follows: Section 2 gives an overview of the MWS data from DESI and introduces the data selections. We introduce the methodology adopted in this paper in Section 3, including the selection function of DESI and forward modeling. Section 4 presents the global best-fit model of the MW stellar halo, the shape and orientation dependence on distances and alignments with the disk and distant satellites, the anisotropic density distribution induced by GSE and LMC, and the model dependence on [Fe/H]. We discuss and compare our results with previous studies in Section 5, and we conclude in Section 6. 2. Data 2.1. The DESI Milky Way Survey In this paper, we use K giants observed in the MAIN- BRIGHT program of DESI MWS DR2 observation, which is the current internal data release and will be released in 2027. The MAIN-BRIGHT program is one of the main target categories of DESI MWS, which is observed in bright time together with DESI Bright Galaxies (C. Hahn et al. 2023) and covers a magnitude range of 16< r < 19. For the MAIN-BRIGHT program of DESI, the targets are further split into a few different classes, including MAIN- BLUE and MAIN-RED. The MAIN-BLUE subsample contains blue stars (g − r < 0.7), which mainly probe metal-poor thick- disk/halo turnoff stars. The MAIN-RED subsample contains red stars (g − r � 0.7), selected with additional astrometric criteria to prioritize the selection and observation of distant halo giants. MAIN-BLUE has the same observation priority as MAIN-RED. In this work, we only use MAIN-BLUE and MAIN-RED subsamples. We refer readers to A. P. Cooper et al. (2023) for more details on DESI MWS target categories and criteria. 2.2. MWS Stellar Parameter Catalog In this work, we use the stellar parameter and radial velocity value-added catalog (VAC) produced by the internal RVJ pipeline, which is an improvement on the DESI MWS RVSpecFit (RVS) pipeline used to produce the early data release (EDR) and the first data release (DR1) VACs (S. E. Koposov et al. 2024, 2025). The biggest difference between the RVJ and RVS pipelines is that the RVJ pipeline uses templates that are synthesized with Korg (A. J. Wheeler et al. 2023) and provides measurements for individual elemental abundances, including [Mg/Fe], rather than just bulk alpha abundances. Distances of individual stars are estimated from the DESI SpecDis2 VAC (Li et al. 2026, in preparation), which is an updated version of the DESI DR1 SpecDis VAC (S. Li et al. 2025). The main difference between the SpecDis2 and SpecDis catalogs is that SpecDis2 predicts distances from stellar spectra and stellar colors, while SpecDis predicts distance modulus only from stellar spectra. Moreover, Spec- Dis2 trains main-sequence and giant stars separately. Spec- Dis2 has significantly increased the precision of distance measurements for giants. The median distance uncertainty of giants decreases from 40% (SpecDis) to 12% (SpecDis2). 2.3. Selection of K Giants In this subsection, we introduce our K giant sample selection criteria. Our goal is to build a highly pure halo K giant sample to allow robust measurements of the shape, orientation, and radial density profile of the stellar halo. Based on the DESI MWS MAIN-BLUE and MAIN-RED subsamples, we first select K-type stars with colors in the range of 0.5< g − r < 1.3 after extinction correction. With this selection, we get a sample that contains ∼135,000 stars. We further select K giants by absolute magnitude (Y. Huang et al. 2023). Our K giants are defined by −3< Mr < 2 (r-band absolute magnitude after extinction correction), which removes various types of faint stars. After this selection, the remaining sample contains ∼92,000 stars.41 Most of these K giants selected through the steps above lie at >1 kpc above the disk plane, which is due to the flux limit of the survey (16< r < 19) and the magnitude cut we adopt here, which has already removed the majority of contamina- tion from the disk. We further remove stars with [Fe/H] > – 0.85 dex (see red dashed line in Figure 1). This may also remove a small number of metal-rich GSE stars (see the lower right corner of Figure 1), but it ensures a purer halo star sample (M. R. Hayden et al. 2015). After this selection, we retain about 46,000 stars. As the final requirement in our selection, we remove stars belonging to known substructures. Member stars from globular clusters (GCs) and dwarf galaxies are removed by matching specific stars with H. Baumgardt & E. Vasiliev (2021) and A. B. Pace et al. (2022). The contribution from the Sagittarius (Sgr) stream is removed by footprint. We compute the latitude coordinate, B̃, perpendicular to the plane of the Sagittarius stream according to the coordinate system defined by V. Belokurov et al. (2013) for all stars in the sample, and we remove stars with B 15˜ < . After this, we get a final halo K giant sample that contains ∼28,000 stars, which covers a galactocentric distance (rGC) from 8 to 200 kpc. 41 We note that MAIN-RED targets are incomplete when close to the Sun, because stars with ω > 3σω +0.3 mas are not included as MAIN-RED targets. However, all K giants in our final sample are beyond 6 kpc from the Sun owing to our chosen absolute magnitude range above. The incompleteness of MAIN-RED targets would not affect our results. 3 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. 3. Methods We take two major steps to construct the spatial distribution of the stellar halo with K giants: (1) correcting the selection effects due to incompleteness of the spectroscopic survey and incomplete- ness of photometric targets due to the flux limit, and (2) re- constructing the stellar halo with forward modeling, which contains a set of flexible and smooth functions with free parameters. 3.1. Correction of Incompleteness Our K giant sample is incomplete for two reasons. The first is incompleteness in spectroscopically observed stars with respect to photometric targets. The second is a consequence of the MWS target selection, as the MAIN-BRIGHT program only observed stars with 16< r < 19, and thus K giants with different ranges of absolute magnitudes are selected at different distances, resulting in an incompleteness that is radially dependent. We correct these two selection effects with two selection functions. 3.1.1. Correcting Incompleteness of the Spectroscopic Survey The first selection function is defined as the completeness fraction of spectroscopically observed stars with respect to the photometric targets, which can be calculated by comparing the spectroscopically observed stars with the photometric input targets. In our final halo K giant sample, stars with 0.5< g − r < 0.7 and 0.7< g − r < 1.3 belong to MAIN- BLUE and MAIN-RED, respectively. We then construct two selection functions for MAIN-BLUE and MAIN-RED sepa- rately.42 According to C. Liu et al. (2017), the selection function (S) usually depends on the sky position, color, and apparent magnitude, while we further check that for both MAIN-BLUE and MAIN-RED the selection function in fact shows almost no dependence on color or apparent magnitude. Thus, for both MAIN-BLUE and MAIN-RED sources, the selection functions reduce to being only dependent on sky coordinates: S n n obs R.A.,decl.,class R.A.,decl. R.A.,decl. . 1 sp,class ph,class ( ) ( ) ( ) ( )= Here “class” refers to MAIN-BLUE or MAIN-RED, “obs” stands for observed by DESI, and R.A., decl. are two celestial coordinates. Parameters nsp and nph are the number of spectroscopically observed stars by DESI and the photometric targets, respectively, as a function of sky coordinates. We call this the angular selection function. To calculate the angular selection function, we first divide the whole footprint of DESI DR2 MWS data into different regions, and each region has a different number of passes, called Npasses. This parameter indicates the number of DESI tiles at a given sky position.43 Naturally, Npasses is a function of R.A. and decl., and a higher Npasses will lead to a higher completeness fraction and thus a higher selection function value. For regions with the same Npasses, we pixelize them with the Python package healpy (A. Zonca et al. 2019) and count the number of spectroscopic and photometric stars in each pixel to get the selection function in different pixels. As an example, we present the selection function of the MAIN-BLUE subsample. Figure 2 shows the angular selection function of MAIN-BLUE as a function of sky position. The tiling patterns are well captured. 3.1.2. Correcting Incompleteness of Photometric Targets In this subsection, we introduce the second selection function that corrects the incompleteness of photometric targets due to the survey flux limit. Since the flux limit affects the distance range that we can observe for a source with fixed absolute magnitude, we call it the radial selection function. We assume that the luminosity function of K giants obeys an exponential law M 10r M0.32 r( ) (X.-X. Xue et al. 2014). For a given distance D(kpc) to the Sun, the apparent magnitude of K giants covers a range from r D7 5 log= + (corresp- onding to an absolute magnitude of Mr = −3) to r D12 5 log= + (corresponding to an absolute magnitude of Mr = 2). The fraction of a K giant that has been targeted is S D m m r D dr r D dr obs , , , class 5 log 10 5 log 10 , 2m m D D radial 1 2 7 5 log 12 5 log 1 2 ( ) ( ) ( ) ( )= + + where m1 is the larger of 16 and D7 5 log+ , and m2 is the smaller of 19 and D12 5 log+ . Again, “class” refers to MAIN-BLUE or MAIN-RED. The typical value of Sradial is about 0.4 at D ∼ 30 kpc. Sradial decreases quickly with distance, which becomes ∼0.03 at D ∼ 50 kpc. With the angular and radial selection functions, we can correct for the incompleteness in our sample of K giants. This is achieved by assigning each star a weight, which is the product of the inverse of the angular and radial selection functions. 2.0 1.5 1.0 0.5 [Fe/H] 0.1 0.0 0.1 0.2 0.3 0.4 0.5 [M g/ Fe ] Halo ThinDisk ThickDisk MetalRich GSE 100 101 102 Co un ts Figure 1. DESI DR2 giants in [Fe/H]–[Mg/Fe] space from the RVJ pipeline. The giants are selected from the MAIN-BRIGHT program. Halo stars are defined by [Fe/H] < –0.85 (vertical red dashed line). 42 It is important to treat MAIN-BLUE and MAIN-RED separately, as they correspond to different photometric targets. 43 For DESI bright time observation, the tiles are organized into five passes. Each tile corresponds to a single arrangement of DESI fibers to observe sources. Tiles within a pass do not overlap, whereas tiles on different passes are offset such that all points on the sky will be covered by three tiles on average (A. P. Cooper et al. 2023), i.e., the mean Npasses is 3. However, most MWS sources are observed only once, as bright galaxies are also observed at bright time and are assigned higher fiber priorities than most MWS sources. 4 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. 3.2. Constructing the Stellar Halo with Forward Modeling In this subsection, we introduce our parameterized forward model, which describes the shape, orientation, and radial density profile of the MW stellar halo. We assume that the stellar halo is a triaxial ellipsoid. The shape of the ellipsoid can be fully described by the flattened radius: r X Y p Z q , 3q 2 2 2 ( )+ + where X, Y, and Z are Cartesian coordinates centered on the galactic center and p, q are intermediate-axis flattening and minor-axis flattening (0< q < p < 1), respectively. Here we allow the X-, Y-, and Z-axes to be rotated with respect to the standard galactocentric Cartesian coordinates (x, y, z), with z perpendicular to the Galactic disk, and x points from the Sun toward the Galactic center. Hence, in principle, we need three rotation angles to fix the orientation of the ellipsoid in our modeling: the rotation angle around the minor axis (f, yaw angle), the rotation angle around the intermediate axis (θ, pitch angle), and the rotation angle around the major axis (ζ, roll angle). In this work, we adopt the convention that a positive value of the rotation angle indicates clockwise rotation. Importantly, the pitch angle θ is the angle between the major axis of the stellar halo and the Galactic plane. The yaw angle f is the angle between the vector connecting the Sun to the Galactic center and the projection of the model major axis on the Galactic plane. However, J. J. Han et al. (2022) found a prolate-like stellar halo and hence failed to get a good fitting of the roll angle. In this work, we also find that the stellar halo is nearly prolate with a similar axis ratio to that in J. J. Han et al. (2022; see Section 4). We also try to fit the roll angle, and the best-fit value is nearly zero. Hence, in this work we choose to only fit the pitch angle and yaw angle of the stellar halo, by fixing the roll angle to zero throughout our analysis. Following J. J. Han et al. (2022), we assume that the radial density profile of the MW stellar halo obeys a triple power law with two break radii r r r r r r r r r r r r r ; , ; , ; , 4q q q b b q b q b b q b q ,1 ,1 ,1 ,2 ,2 ,2 1 2 1 2 3 2 3 ( ) ( )< < where α1, α2, and α3 are the first, second, and third power-law slopes, respectively, and rb,1 and rb,2 are the first and second break radii, respectively. We also test a less flexible radial density profile with only one break radius and two power-law slopes, r r r r r r r r ; , ; , 5q q q b b q b q ,1 ,1 ,1 1 2 1 2( ) ( )< to compare with the triple power-law model. In total, there are nine free parameters for the triple power-law model (two in shape, two in orientation, and five in density profile) and seven for the second power-law model, which lacks two density profile parameters. We have also tried other, more asympto- tically changing functional forms, but we have failed to see significant changes in the best-recovered positions of the break radii. In addition to the above triple power-law model, we also adopt a model that allows the shape and orientation to change with distances, which is called the variable model.44 We assume that the p, q, θ, and f parameters change with galactocentric radius, rGC, as third-order polynomials: p r p k r k r k r , 6p p pGC 0 1 GC 2 GC 2 3 GC 3( ) ( )= + + + q r q k r k r k r , 7q q qGC 0 1 GC 2 GC 2 3 GC 3( ) ( )= + + + r k r k r k r , 8GC 0 1 GC 2 GC 2 3 GC 3( ) ( )= + + + r k r k r k r . 9GC 0 1 GC 2 GC 2 3 GC 3( ) ( )= + + + Here p0, q0, f0, and θ0 are axis flattening parameters and rotation angles at rGC = 0. Parameters kpi, kqi, kfi, and kθi are 0 50 100 150 200 250 300 350 ra 40 20 0 20 40 60 80 de c Excluded Sgr region 0.0 0.2 0.4 0.6 0.8 1.0 S a ng ul ar Figure 2. Angular selection function of the MAIN-BLUE subsample as a function of sky position. The sky region between the two red dashed lines denotes the Sgr footprint, which has been cut off from our analysis. Here the angular selection function is defined as the completeness fraction of DESI-observed stars with respect to the targets. Brighter colors indicate higher completeness. The MAIN-RED subsample has a similar angular selection function to MAIN-BLUE. We do not repeatedly present the angular selection function of the MAIN-RED subsample. 44 For the variable model, we fix the radial density profile to be a triple power- law model, with parameters from those of the best-fit model obtained in Section 4.1. 5 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. polynomial coefficients, and i is in the range of 1−3. There are 16 free parameters in the variable model. We will present the best-fit triple power-law model in Section 4.1 and compare it with the best-fit double power-law model. We further show the best-fit variable model in Section 4.2. 3.3. Likelihood Function With angular and radial selection functions and a forward model for the radial density profile above, the probability for the ith K giant with a given heliocentric distance Di and sky coordinates R.A.i and decl.i to be observed by DESI is P D S S r r C , R.A. , decl. 4 , 10i i i i i i q i q i,angular ,radial , , 2 ( ) ( ) ( )= C S S r D d dD, 11 D qangular radial 2( ) ( )= where C is the normalization and integrates over the volume in the (D, R.A., decl.) space. Si,angular and Si,radial refer to the selection functions defined in Equations (1) and (2), respec- tively. Equation (10) takes the same form for stars in either MAIN_BLUE or MAIN_RED, adopting the corresponding MAIN_BLUE or MAIN_RED selection functions, respectively. Then, we can directly calculate the likelihood of observed stars L P DR.A. , decl. , , 12 i N i i i i 1 K ( ) ( )= = where NK is the total number of K giants and thus the multiplication goes through all K giants in both MAIN_BLUE and MAIN_RED. We then assume a uniform prior for all parameters in the stellar halo model and use the Python package emcee (D. Foreman-Mackey et al. 2013) to sample the posterior distribution. We employ 50 walkers, discard the first 400 steps as before burn-in,45 and run 1000 steps in total. 4. Results In this section, we first present the results of the best-fit triple power-law model for the MW stellar halo (see Section 4.1). We then show the results of the best-fit variable model (see Section 4.2). We further investigate the spatially anisotropic stellar halo by measuring the density profiles along different pointings on the sky, including regions associated with the GSE debris and the LMC wake (see Section 4.3). Finally, we explore the model dependence on [Fe/H] by fitting the triple power-law model in various [Fe/H] ranges (see Section 4.4). We have verified in the Appendix that distance uncertainties of K giants do not significantly impact our results. 4.1. Best-fit Triple Power-law Model Figure 3 shows the posterior contours of each parameter of the triple power-law model. The best-fit parameters and their associated uncertainties are shown in the first row of Table 1. For our best-fit model, we find the following: (1) The axis ratios of the stellar halo are about 1:p:q = 10:8:7. Hence, the stellar halo is nearly prolate with p and q more similar to each other. (2) The stellar halo is tilted by about 44° off the Galactic plane and about 27° away from the Sun−Galactic center axis. (3) The stellar halo has two break radii at 16 and 76 kpc (see inset panel of Figure 3). Figure 4 shows the comparison between our best-fit parameters of the triple power-law model (Equation (4)) and the best fits reported by J. J. Han et al. (2022). Our shape parameters p and q are similar to those of J. J. Han et al. (2022). We also found a similar yaw angle f, while the pitch angle θ is more negative than that of J. J. Han et al. (2022). The discrepancy in orientation parameters may be explained as follows: (1) DESI and H3 surveys have different footprints. (2) J. J. Han et al. (2022) removed Sgr by angular momenta,46 while we mask the footprint of the Sgr stream, which adds a difference in footprint. However, we find a more prominent discrepancy in the radial density profile. The first break radius rb,1 is larger than that in J. J. Han et al. (2022), and we do not find a break radius around 30 kpc like J. J. Han et al. (2022), but instead we find a much greater second break radius rb,2 at ∼70 kpc. Our radial density profile is also steeper than that of J. J. Han et al. (2022) beyond the first break radius. The discrepancy in the radial density profile is mainly caused by the different data selection compared with J. J. Han et al. (2022). J. J. Han et al. (2022) mainly focused on the radial density profile of GSE stars, by adopting chemical and kinematic selections to obtain a pure GSE sample. Their sample is shallower (rGC < 60 kpc), and thus J. J. Han et al. (2022) did not report a break radius at ∼70 kpc. In this work, we use all halo stars. The selection on [Fe/H] (see Section 2.3) further removes some GSE stars, leading to a slightly smaller fraction of GSE stars in our halo sample. However, if we only fit the radial density profile of HAC-N and HAC-S, which are dominated by GSE stars (I. T. Simion et al. 2019; H. D. Perottoni et al. 2022), we can also find a break radius at ∼30 kpc (see Section 4.3). Interestingly, M. Cavieres et al. (2025) utilized a K giant and MSTO sample that can extend to about 100 kpc and also found a break radius at ∼70 kpc in the southern sky. Figure 4 also shows the comparison between the parameters and uncertainties of the best-fit double power-law model (Equation (5)) and triple power-law model. The best-fit double power-law model is also obtained from Markov Chain Monte Carlo (MCMC), for which we do not show the posterior contours. The double power-law model shows similar shape, orientation, and radial density profile parameters except for the second power-law slope α2. The difference between α2 of the two models is greater than the model uncertainties, meaning that a third component is required to model the MW stellar halo. We have also employed the Akaike information criterion47 (AIC; H. Akaike 1974), and the AIC difference between the triple power-law model and the double power-law 45 The autocorrelation time is about 35 steps, so the post-burn-in samples can be treated as statistically independent draws from the target distribution (D. Foreman-Mackey et al. 2013). 46 We have also tested removing Sgr by angular momentum. However, there is still some contamination after the angular momentum selection. Besides, Kizhuprakkat et al. (2026, in preparation) utilize a K giant sample from DESI MWS DR2 for studying the assembly history of the stellar halo. They further confirm that distance errors scatter Sgr members into retrograde and low- energy regions, blending them with GSE, leading to less robust removal of Sgr only by angular momenta. 47 k LAIC 2 2 ln ˆ= , where k is the number of free parameters and Lln ˆ is the maximum likelihood; lower AIC favors a better trade-off between fit quality and model complexity. 6 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. model is −216, which means that the data prefer the triple power-law model. Figure 5 presents the direct one-dimensional comparison between the observed and the best-fit double/triple power-law model profiles, reported as a function of the flattened radius, rq. As a comparison, we also plot the best-fit profile from J. J. Han et al. (2022). In the top panel, the green circles present the observed stellar density48 profile multiplied by rq 2. The error bars give 1σ uncertainties estimated from 100 bootstrap resamplings of our K giants. The red and blue lines represent the predicted density profile from the triple and double power- law models, for which we convolve them with the DESI footprint and selection functions. The bottom panel of Figure 5 presents the ratio between the best-fit model profiles and the real data as a function of rq. Both the triple and double power-law models represent good fits to the observation within 70 kpc. The differences between the two models and the observed density profile are all smaller than 10%, with the double power-law model showing slightly larger bias. However, the double power-law model fails to fit the observed density profile and predicts an overestimated density profile beyond 70 kpc, while the triple power-law model captures the transition of the density profile with a second break radius rb,2 around 70 kpc. The triple power-law model better describes the radial density profile of the stellar halo than the double power-law model beyond 70 kpc. 4.2. Radially Variable Model with Fixed Radial Density Profile In this subsection, we explore a more flexible radially variable stellar halo model, which allows the shape and orientation parameters to change with the radius (see details in Section 3.2). We only explore this model with the whole sample, as Figure 12 above has shown that the dependence of shape and orientation parameters on [Fe/H] is relatively weak. We adopt the radial density profile to be a fixed triple power-law model, with the parameters rb,1, rb,2, α1, α2, α3 fixed to best-fit values from Section 4.1, and then introduce the dependence of the shape and orientation parameters on radius (see Equation (6)). Figure 6 shows the posterior contours for the variable model. Different orders of kfi/kθi are highly correlated with each other, whereas the other model parameters are more independently constrained. In the upper right corner of Figure 6, we present two panels that show the dependence of Figure 3. Posterior contours for different combinations of nine parameters of the triple power-law model. Red and blue dashed lines indicate 50th, 16th, and 84th percentiles. The meanings of different model parameters can be found in Section 3.2. The yellow, light-blue, and dark-blue contours represent the 30%, 1σ, and 2σ regions of the MCMC post-burn distributions, respectively. The black line in the inset panel shows the selection-effect-free best-fit model radial density profile (renormalized to unity), with the red shaded region representing the 1σ model uncertainty, and two vertical black dashed lines mark two break radii. 48 The density is obtained from dividing the number of stars in each ellipsoidal shell by its volume. 7 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. Table 1 Best-fit Parameters of Triple Power-law Model in Different Regions Name (l, b) rb,1 rb,2 α1 α2 α3 p q θ f (kpc) (kpc) (deg) (deg) Halo 0° < l < 360°, −90° < b < 90° 16.0 0.20 0.25+ 76.3 1.94 2.37+ 1.50 0.07 0.07+ 3.45 0.02 0.02+ 5.20 0.14 0.14+ 0.85 0.01 0.01+ 0.74 0.01 0.01+ 43.8 0.7 0.7+ 26.9 1.0 1.0+ P1 0° < l < 180°, 0° < b < 90° 14.4 0.2 0.2+ 69.3 3.5 3.3+ 1.28 0.08 0.07+ 3.48 0.02 0.02+ 4.84 0.16 0.17+ 0.93 0.01 0.01+ 0.71 0.01 0.01+ 38.6 0.8 0.8+ 0.4 0.6 0.3+ P2 180° < l < 360°, 0° < b < 90° 17.7 1.0 0.7+ 103.1 2.2 2.9+ 1.72 0.33 0.17+ 3.56 0.02 0.03+ 7.47 0.15 0.47+ 0.71 0.07 0.08+ 0.66 0.03 0.03+ 52.9 1.0 1.0+ 14.2 8.0 8.0+ P3 0° < l < 180°, −90° < b < 0° 14.9 0.5 0.5+ 60.0 0.9 3.0+ 1.34 0.15 0.15+ 3.25 0.03 0.03+ 4.89 0.15 0.15+ 0.66 0.01 0.01+ 0.64 0.01 0.01+ 97.0 0.7 0.7+ 30.3 1.6 1.6+ P4 180° < l < 360°, −90° < b < 0° 17.9 10.6 1.5+ ⋯ 0.91 0.39 1.57+ 3.89 1.51 0.16+ ⋯ 0.77 0.08 0.04+ 0.73 0.03 0.06+ 5.8 4.5 3.4+ 1.6 2.6 1.3+ Pisces (l − 74°)2 + (b + 47°)2 < 400 9.5 4.5 5.1+ 60.1 10.1 5.3+ 1.95 0.99 1.10+ 3.22 0.04 0.06+ 4.99 0.18 0.25+ 0.74 0.03 0.04+ 0.72 0.03 0.06+ 76.1 3.1 3.8+ 21.6 63.6 10.0+ HAC-S 30° < l < 60°, − 45° < b < − 20° 34.3 8.1 5.8+ ⋯ 0.52 0.48 0.31+ 2.63 0.03 0.05+ ⋯ 0.21 0.02 0.06+ 0.17 0.03 0.05+ 100.1 1.2 1.4+ 0.3 0.44 0.3+ HAC-N 30° < l < 90°, 20° < b < 45° 31.8 6.4 9.5+ ⋯ 1.71 0.15 0.09+ 2.98 0.06 0.06+ ⋯ 0.32 0.06 0.09+ 0.31 0.06 0.09+ 52.8 2.7 2.8+ 0.9 1.1 0.6+ VOD 270° < l < 330°, 50° < b < 75° 15.8 2.6 8.6+ ⋯ 1.42 1.06 0.98+ 3.15 0.22 0.27+ ⋯ 0.69 0.21 0.18+ 0.59 0.09 0.10+ 55.1 9.9 13.8+ 10.1 20.5 8.2+ 8 T h e A stroph ysical Jou rn al, 999:108 (19pp), 2026 M arch 1 L i et al. rotation angles (top panel) and flattenings (bottom panel) on rGC. The thick lines are the predictions from the best-fit model, and the shaded regions present the 1σ uncertainty of the variable model. We find that both the pitch and yaw angles, θ and f, decrease with increasing rGC. The yaw angle is more aligned with the line connecting the Sun and Galactic center within 20 kpc, which increases to larger absolute values at larger distances. Interestingly, there is a clear trend that the inner stellar halo within rGC≲ 30 kpc is more aligned with the disk, with the absolute value of the pitch angle, θ, smaller than 45°. However, the absolute value of θ becomes greater than 45° beyond 30 kpc, indicating that the major axis of the outer stellar halo becomes more perpendicular to the MW disk, i.e., the orientation of the stellar halo flips with distance. Here a negative pitch angle means that the major axis of the stellar halo is pointing to the south Galactic pole. The pitch angle even reaches −80° when rGC is around 60 kpc. Moreover, the flattening parameter, q, of the minor axis is first increasing with the radius, reaching a maximum at rGC ∼ 45 kpc and then decreasing again. The flattening of the intermediate axis, p, keeps decreasing. In fact, the stellar halo is more prolate at rGC > 30 kpc, where p and q are closer to each other, while both are smaller than 1. At a smaller rGC, p is closer to 1, whereas q is smaller by ∼0.15, indicating that the more nearby inner stellar halo is more oblate. Figure 7 more intuitively shows how the best-fit variable ellipsoidal model changes with radius in galactocentric Cartesian coordinates. The red solid lines are isodensity contours. The black dashed lines indicate reference circles. The left and middle panels show the density distributions in the x–z (y = 0) and y–z (x = 0) slices that are perpendicular to the disk. The stellar halo exhibits a pronounced twist: the isodensity contours systematically rotate toward perpendicular to the disk plane, and the flattening changes with increasing galactocentric distance. The right panel shows the density distribution in the disk plane (z = 0). We also see a trend that on this plane the stellar halo is closer to circles at smaller radii, which becomes more elliptical with the increase in radius. Combining the trends in orientation angles and the two axis ratios with rGC, it is straightforward to see that the inner stellar halo within 30 kpc is oblate and the major axis is more aligned with the disk (with the minor axis perpendicular to the disk). On the other hand, the more distant outer stellar halo becomes more prolate and perpendicular to the disk, with the major axis more aligned with the VPOS of MW satellites, GCs, and stellar streams (e.g., M. S. Pawlowski et al. 2012; M. S. Pawlowski & P. Kroupa 2014). However, the trends of flattening parameter and orientation angles at a larger distance (rGC > ∼ 50–60 kpc) are still very uncertain owing to quite large model uncertainties. Similar to Figure 5, Figure 8 presents the comparison between the observation, the best-fit triple power-law model, and the variable model profile, as a function of the flattened radius, rq. The two vertical black dashed lines present two break radii from the best-fit triple power-law model in the previous subsection. Both the triple power-law model and variable model fit the one-dimensional density profile well at rb, 1 12 14 16 rb, 2 30 40 50 60 70 80 1 1.4 1.6 1.8 2 3.0 3.2 3.4 3.6 3 4.6 4.8 5.0 5.2 p 0.78 0.80 0.82 0.84 0.86 q 0.71 0.72 0.73 0.74 0.75 45 40 35 30 25 30 28 26 24 22 20 Double power-law Triple power-law Han et al. 2022 Figure 4. Best-fit model parameters and uncertainties for the double power-law model (red), the triple power-law model in this work (black), and the triple power- law model in J. J. Han et al. (2022; green). The error bars represent the 1σ uncertainties of three models. The meanings of different model parameters can be found in Section 3.2. Three models share similar flattening parameters and yaw angle (the difference is within 1σ uncertainty), while they show a greater discrepancy in radial density profile parameters. 106 107 108 r2 q (r q )[k pc 1 ] Triple power-law Double power-law Han_2022 Observation 101 102 rq[kpc] 0.5 1.0 1.5 M od el /D at a Figure 5. The top panel shows with green circles the observed stellar density multiplied by rq 2 as a function of flattened radius, rq, compared to the best-fit triple power-law model (red solid line) and double power-law model (blue solid line). Error bars represent the 1σ uncertainties computed from 100 bootstrap subsamples of our K giants. Here the best-fit triple power-law model in this work (red solid line), that in J. J. Han et al. (2022; orange solid line), and the double power-law model (blue solid line) have been convolved with the angular and radial selection functions to have a fair direct comparison with the data (see Section 3.1). Two vertical black dotted lines in both the top and bottom panels represent the two break radii, rb,1 and rb,2, respectively. Two horizontal dashed gray lines in the bottom panels represent 10% regions of the model-predicted density over observed density, with the black dashed horizontal line marking y = 1. Both triple power-law and double power-law models fit the stellar halo within 70 kpc well, but the double power-law model shows a worse match of the outer stellar halo beyond ∼70 kpc. 9 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. rb,1 < rq < rb,2, while the variable model fits the observation better within the first break radius. It is hard to distinguish which model fits the observation better at rq > rb,2, due to the much smaller number of stars there and the large error bars. Nevertheless, the difference of the AIC values between the variable model and the triple power-law model is about −136, which means that the data prefer the variable model. We have further checked the angular residuals of the two models and found that the variable model has smaller residuals. Figure 9 shows the angular residuals of the densities predicted by the variable model in Mollweide projection in Galactic coordinates and in four different ranges of rGC (see the text on top of each panel). The density at a given pixel is the ratio between the number of stars at this pixel and the corresponding volume. The angular residual is further defined by the difference between observed density and predicted density at different l and b, divided by the prediction at a given distance range. Here the LMC is marked by a black pentagram. However, we find that, even with the most flexible variable model, we still cannot fully capture all angular substructures of the real stellar halo. There are some overdense regions in Figure 9. For example, we find that the region enclosed by the red solid line presents some overdensity in all four panels at different distances. We also observe a mild overdensity of the areas enclosed by the green, blue, and yellow solid lines, although it is not present in every panel. The enclosed regions by the red, yellow, and blue solid lines are associated with three well-known overdense Figure 6. The posterior contours of 16 model parameters of the variable model. Red dashed lines and blue dashed lines indicate 50th, 16th, and 84th percentiles. The meanings of different model parameters can be found in Section 3.2. The yellow, light-blue, and dark-blue contours represent the 30%, 1σ, and 2σ regions of the MCMC post-burn distributions, respectively. The two panels in the upper right corner present the dependence of the rotation angles f, θ and flattenings p, q on rGC. The thick lines are the predictions from the best-fit model, and the shaded regions present the 1σ uncertainty of the variable model. 10 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. regions, HAC-N, HAC-S, and VOD (A. K. Vivas et al. 2001; H. J. Newberg et al. 2002; V. Belokurov et al. 2007; M. Jurić et al. 2008; L. L. Watkins et al. 2009; B. Sesar et al. 2010; I. T. Simion et al. 2014), which are chemodynamically consistent with GSE stars (I. T. Simion et al. 2019; R. P. Naidu et al. 2021; H. D. Perottoni et al. 2022). Besides, the region enclosed by the green line contains Pisces overdensity (V. Belokurov et al. 2019; C. Conroy et al. 2021), which could be associated with the transient density wake of the LMC. The definition of these solid-line-enclosed regions can be found in Table 1. We also see some slightly underdense regions in each panel. These slightly underdense regions might be caused by the existence of the overdense regions we talked about before, and thus the model tends to be shifted to larger values in order to better match the few overdense regions. We further explore these substructures in Section 4.3. 4.3. Spatial Anisotropies of Stellar Halo Density Distribution and the LMC Density Wake In Sections 4.1 and 4.2, we fit the stellar halo with a continuous and smooth ellipsoid model. However, our MW stellar halo is full of substructures, and even the variable model cannot fully capture these density fluctuations. Thus, we show the observed density contrast in Figure 10 in the Mollweide projection in Galactic coordinates. The density contrast is defined by the difference between the observed density at given (l, b) and the mean observed density in the ellipsoidal shell at the given distance and is divided by the latter. To get rid of any artificial over- or underdense features, we have corrected the incompleteness with angular and radial selection functions (Equations (1) and (2)). There are a few overdense regions, marked by colored solid lines, which are roughly revealed in Figure 9. In addition to these overdense regions, there is an overdense region in the upper right part of the sky (210° < l < 330° and b > 0°), with the overdensity becoming more prominent with the increase in distances. For these overdense regions, the stellar halo may have different break radii and power-law slopes. We then divide the footprint into a few pieces and model the shape, orientation, and radial density profile separately. We explore the aniso- tropy in eight regions defined in Table 1 with the triple power- law model. Here HAC-N,49 HAC-S, and VOD are spatial overdensities observed in the sky and chemodynamically associated with the GSE debris (I. T. Simion et al. 2019; H. D. Perottoni et al. 2022). Pisces overdensity is correlated with the LMC50 (N. Garavito-Camargo et al. 2019; C. Conroy et al. 2021; M. Cavieres et al. 2025). P1, P2, P3, and P4 are the first, second, third, and fourth quadrants in the galactic projection. P2 is an overdense region, while P4 is an underdense region. These are likely associated with the global 50 0 50 x[kpc] -50 0 50 z[ kp c] 50 0 50 y[kpc] -50 0 50 z[ kp c] 50 0 50 x[kpc] -50 0 50 y[ kp c] Figure 7. The density distribution of the best-fit variable model in the x–z (y = 0), y–z (x = 0), and x–y (z = 0) slices. Here x, y, z are standard Galactocentric Cartesian coordinates. The red solid lines are isodensity contours, with the density for each outer contour decreased by a factor of 10 compared with the previous inner one. The black dashed lines indicate reference circles. The yellow pentagram denotes the location of the Sun. 106 107 108 r2 q (r q )[k pc 1 ] Triple power-law Variable model Observation 101 102 rq[kpc] 0.5 1.0 1.5 M od el /D at a Figure 8. Similar to Figure 5, but the blue solid curve shows the variable model, by allowing the orientation angles and axis ratios to vary with radius. The red solid curve and green circles with error bars are the same as those in Figure 5, which are the best-fit triple power-law model profile and the real data. The model predictions have been convolved with the angular and radial selection functions (see Section 3.1) and have been renormalized to give the same total number of stars as real data. 49 Note that the area of HAC-N defined in this work is twice as large as previous studies (J. A. S. Amarante et al. 2024). We have further checked that the area (30° < l < 60°, −45° < b < −20°) and HAC-N defined in this work share similar radial density profile, shape, and orientation parameters by forward modeling. Hence, we believe that a different definition of HAC-N would not affect the results and keep the definition of HAC-N in the rest of this work. 50 N. Garavito-Camargo et al. (2019) presented high-resolution N-body simulations of the MW−LMC interaction. They identified two main density wakes, the transient wake and the collective wake in the MW stellar halo induced by the LMC. The transient wake is a trailing overdensity following the orbit of the LMC, while the collective wake is a broader, persistent overdensity in the northern halo caused by resonant global response. 11 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. collective response of the stellar halo to LMC infall, and we will discuss them later in this section. Figure 11 presents the comparison between the global best- fit variable model obtained in Section 4.2 and the corresp- onding triple power-law models (different color curves in each panel) that directly fit the measured density profile in each region (see the legend). The global variable model profiles differ in different panels because the variable model is dependent on angular position. Here we have corrected the selection effect. In all panels of Figure 11, the triple power-law model at a given sky position fits the observed density profile well. We also find that the variable model does not present a very large difference compared with the best-fit triple power- law model in each region, which indicates that the variable model is flexible in terms of capturing the main density variations over the sky. However, in the HAC-N region, observation points are slightly higher than the global black curve at all distances, while in the HAC-S and Pisces regions the observation points are slightly higher than the black curve at smaller distances. And in the P4 region, the variable model slightly overestimates the density profile. These are all consistent with Figure 9. Table 1 provides the best-fit parameters in each panel of Figure 11. In particular, we find that for the best-fit triple power-law model of HAC-S, HAC-N, VOD, and P4 regions, the difference between the best-fit second power-law slope α2 and the third power-law slope α3 is negligible when compared with the uncertainties of α2 and α3. Hence, the triple power- law model degrades to a double power-law model, and thus we omit the second break radius rb,2 and α3 in these regions. The break radius of VOD is about 15 kpc, which is consistent with the first break radius we find in Section 4.1. The break radii of both HAC-N and HAC-S are around 30 kpc, which is consistent with the second break radius reported in J. J. Han et al. (2022) for GSE. This may indicate that HAC-N and HAC-S are highly correlated with GSE stars (I. T. Simion et al. 2019; H. D. Perottoni et al. 2022; D. Ye et al. 2024). Besides, we find that the ellipsoids of HAC-N and HAC-S are highly prolate and present very large pitch angles. We report a pitch angle of 100° and 53° for HAC-S and HAC-N, respectively. However, given the small area of these overdense regions, the best-constrained axis ratios, pitch angles, and yaw angles could be subject to large uncertainties due to extrapolations, though the associated statistical errors are small, and thus we do not overinterpret these shape and position parameters. Moreover, we have tried directly fitting spherical models. The best-fit density slope parameters and break radii remain consistent with those of the ellipsoidal models adopted here. For P1 and P3 regions, which contain HAC-N and HAC-S, respectively, the best-fit radial density profile is a triple power- law model. P3 is nearly prolate, while P1 is more likely to be oblate. P1 and P3 share similar break radii, with rb,1 around 14 kpc and rb,2 around 60−70 kpc, which is consistent with Section 4.1. We also see a second break radius around 60 kpc in the Pisces region, which is consistent with M. Cavieres et al. (2025), who measured the radial density of the Pisces region based on MSTO stars and found a break radius at ∼70 kpc. P3 has a similar second break radius to Pisces, which is likely because the P3 region is dominated by Pisces at larger distances. N-body models (e.g., N. Garavito-Camargo et al. 2019; C. Conroy et al. 2021; Y. Sheng et al. 2024, 2025; V. Chandra et al. 2025) predict that the first infall of the LMC produces an overdensity at ∼60−70 kpc in the Pisces region, and thus the second break radius we found in Pisces might be related to the LMC transient wake. For the P2 and P4 regions, the best-fit density profiles deviate from those in the P1 and P3 regions. We find that the second break radius of P2 is approximately 100 kpc. More- over, we see a clear overdensity of the observational dots, peaked at 90 kpc, as marked by the red triangular box with the zoom-in figure shown in the lower left corner of the P2 panel. The obvious overdensity and a break radius at 90–100 kpc may be new evidences of the LMC collective wake. Moreover, the 306090120150 210240270300330360 LMC VODHAC-N HAC-S Pisces P1 P2 P3 P4 10 – 1.4 dex, but it shows some more prominent change for the last two metal-poor bins. Figure 13 shows a direct comparison among the best-fit model density profiles for subsamples with different metalli- cities. The model density profiles have been renormalized to unity. It is clear that with the decrease in [Fe/H] the distribution of halo stars becomes more extended. Since the break radius is likely correlated with an accretion event (e.g., L. L. Watkins et al. 2009; A. J. Deason et al. 2013; N. Hernitschek et al. 2018; R. P. Naidu et al. 2021; J. J. Han et al. 2022), the dependence of rb,1/rb,2 may reflect a time- dependent stripping process of a dwarf galaxy to form the outer and inner stellar halos at different times. In general, it is expected that the more metal-poor outskirts of a dwarf galaxy are stripped earlier, at larger distances from the center of the host, in contrast to the more metal-rich inner regions, which are later disrupted/stripped at closer distances (H. H. Koppelman et al. 2020; J. A. S. Amarante et al. 2022). Thus, the more metal-poor stellar halo tends to be more extended. Based on Figures 12 and 13, we conclude that the radial density profile has a clear dependence on metallicity, showing that the metal-poor stellar halo is more extended, whereas the Figure 11. In each panel, we compare the predicted density profiles of the variable model from Section 4.2 in the corresponding sky position (black curve) with the result of the triple power-law model (different color curves; see the legend), which directly fits this region. The dashed lines stand for the break radii of the density profile in this region. Note that although the black curves are based on the global fitting with the full data, they differ across different panels, as each panel corresponds to a particular pointing. Selection effects have been corrected in this figure. All colored model density profiles have been renormalized to give the same number of stars as real data at the given sky position after correcting selection effects. The black curve is renormalized to the whole sample. 14 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. dependence of shape and orientation parameters on metallicity also exists but is weaker. 5. Discussions 5.1. Physical Interpretations on the Twisted Inner and Outer Stellar Halos Our results reveal that the inner (≲30 kpc) and outer stellar halos of our MW have different shapes and orientations. The inner stellar halo is oblate and more aligned with the disk, whereas the outer stellar halo is more prolate and perpend- icular to the disk. Our finding here is likely related to the VPOS of our MW (M. S. Pawlowski et al. 2012; M. S. Pawlowski & P. Kroupa 2014), that about 11 classical MW satellites are moving on orbits perpendicular to the MW disk. Nowadays, a growing number of studies argue that groups of satellites moving coherently do exist in ΛCDM simulations, which may form from the same group of galaxies (e.g., M. Cautun et al. 2015a, 2015b; S. Shao et al. 2018, 2019; T. Sawala et al. 2023). Nevertheless, the VPOS indeed shows a group of satellites having high inclination orbits perpendicular to the disk. Now our work here also shows that the outer stellar halo is more prolate and perpendicular to the disk. Hence, the orientation of MW satellites and the outer stellar halo is more consistent, and the MW disk plus inner halo seems to be misaligned with the matter distribution in the outskirts. Thus, the key is to explain the misalignment between the inner and outer matter distribution of our MW. First of all, in numerical simulations, the existence of baryons is predicted to affect the angular momentum of the inner dark matter halo (R. Emami et al. 2021; K. T. E. Chua et al. 2022). P. Bett et al. (2010) pointed out that the formation of the galaxy spins up the dark matter within 0.1 times the virial radius, such that the specific halo angular momentum increases by ∼50% in the median. This results in a better alignment between the angular momenta of the central galaxy and the inner host dark halo. Though we measure the shape and orientation of the MW stellar halo, we may expect to see similar trends between the stellar and dark halos. Another possible scenario could be related to the change or flip in the MW disk and inner halo angular momentum. It was shown by F. A. Gómez et al. (2017) through cosmological hydrodynamical simulation that after the infall of a satellite galaxy the angular momentum of the satellite and the host 14 16 18 r b ,1 [k pc ] All Metallicity error 75 80 85 r b ,2 [k pc ] 0.5 1.0 1.5 1 3.2 3.4 3.6 2 5.0 5.5 6.0 3 0.82 0.84 0.86 0.88 p 1.751.501.251.00 [Fe/H] 0.73 0.74 0.75 0.76 0.77 q 1.751.501.251.00 [Fe/H] 44 42 1.751.501.251.00 [Fe/H] 30 25 20 15 Figure 12. Dependence of the triple power-law model parameters on [Fe/H]. Each green circle presents a subsample of K giants in a given [Fe/H] range, with the error bars showing the 1σ uncertainties of model parameters. The x-axis value of each green circle is the median [Fe/H]. The black horizontal solid line is the best-fit model parameter of the triple power-law model in Section 4.1, with the gray shaded region representing the model uncertainties. 101 102 rq[kpc] 10 7 10 6 10 5 10 4 10 3 10 2 De ns ity [Fe/H]=-0.98 [Fe/H]=-1.05 [Fe/H]=-1.14 [Fe/H]=-1.22 [Fe/H]=-1.34 [Fe/H]=-1.49 [Fe/H]=-1.71 Figure 13. The best-fit selection-effect-free radial density profiles for subsamples with different metallicities (see the legend). All profiles have been renormalized to unity. 15 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. galaxy disk would align. This is not only because the orbital plane of the infalling satellite is affected by the disk potential but also because the disk is responding to the infalling massive satellite. Moreover, studies based on numerical simulations have reported that the spin vector of the inner halo experiences much more frequent flips than the halo as a whole (P. E. Bett & C. S. Frenk 2012), and a significant fraction of halos having had a large spin flip are associated with minor mergers (P. E. Bett & C. S. Frenk 2016). If the above scenario is also true for our MW, it indicates that our MW disk plane and the inner stellar halo may have flipped by ∼90°, due to the infall of a massive satellite. Note that the more frequent flips in the central galaxy and inner halo are likely due to their shorter dynamical timescales than the outer halos. Observational constraints on the MW dark halo shape often rely on stellar streams or halo stars spanning wide ranges of sky footprint. Through the modeling of stellar streams (J. Nibauer & A. Bonaca 2025; C. G. Palau et al. 2025), C. G. Palau & J. Miralda-Escudé (2023) also identified slightly prolate MW dark halos, with their major axis more perpend- icular to the disk. Very recently, starting from the assumption that, given a correct potential model, the distribution function of MW halo stars does not evolve with time (e.g., J. Han et al. 2016a, 2016b; Z. Li et al. 2025), and with the model extended to the nonspherical case (L. Zhu et al. 2025), L. Zhu et al. (2026) reported that the MW dark halo within 50 kpc is nearly oblate, with the long-axis−intermediate-axis plane of the dark matter halo perpendicular to the Galactic disk beyond 20 kpc. Similar conclusions are reported for simulated MW-like systems (e.g., S. Shao et al. 2021). All these studies are showing a broadly consistent picture of our measured MW inner and outer stellar halo orientations. In the end, we would like to emphasize that there is a large scatter in the alignment angles between the spins or major axes of galaxies and host dark halos, with a median misalignment angle of about 30° as reported in both simulations (P. Bett et al. 2010) and real observations by calculating the galaxy shape−shape correlations and shape−shear correlations (e.g., T. Okumura et al. 2009; T. Okumura & Y. P. Jing 2009; K. Xu et al. 2023a, 2023b). Moreover, the correlation strength is much weaker in disk galaxies than in elliptical galaxies (e.g., D. Kirk et al. 2015; K. Xu et al. 2023a). This is likely because galaxy disks are mainly formed from the accreted gas during a certain period of time, whereas the spin of the outer halo is formed by later accreted material, which does not necessarily have the same spin as early accreted gas. 5.2. Comparison with Other Studies Finally, we present the comparison of the radial density profile between this work and the literature and give a brief discussion on the tension in radial density profiles. Figure 14 compares the radial density profiles from this work (red solid line) to those in the literature. The vertical lines represent the best-fit break radii. The horizontal lines stand for the power indices of the radial profiles, with the shaded range presenting the uncertainties. Many works adopt the double power-law model, and the break radii vary from 10 to 30 kpc. We adopt the interpretation from R. P. Naidu et al. (2021) and J. J. Han et al. (2022) and interpret the smaller break radii to be correlated with the apocenter of closer passage of GSE and the greater ones to be associated with the apocenter of the farther passage of GSE. The disparities in measured break radii could be explained by the fact that different surveys have varying sky footprints, as we have shown in Section 4.3 Table 1, the break radii vary with sky positions. We suggest that the break radius beyond 60 kpc is associated with the infall of the LMC (see also Section 4.3). Most previous works failed to find such a great break radius, which is mainly due to the flux limit. For example, the limiting magnitude of the LAMOST survey in Gaia G band is approximately 18 mag (A. L. Luo et al. 2015). C. Yang et al. (2022) only mapped the density distribution of the stellar halo within ∼50 kpc with the LAMOST K giant. And for the H3 survey, the limiting magnitude in Gaia G band is also approximately 18 mag (C. Conroy et al. 2019), and thus J. J. Han et al. (2022) did not report a break radius at 60–70 kpc. There can also be large differences in power-law indices, which are likely due to variations in sky positions and the presence of substructures (e.g., B. Lowing et al. 2015). Table 1 shows that the first power indices can span from 0.5 to 2, the second power indices vary from 2.6 to 3.9, and the third power indices range from 4.8 to 7.5, reflecting significant anisotropies of the stellar halo. The anisotropic stellar halo can be explained by the remnants of major merger events in the sky. N-body simulations of the GSE merger (R. P. Naidu et al. 2021) and the ongoing infall of the LMC (N. Garavito-Camargo et al. 2019; C. Conroy et al. 2021; Y. Sheng et al. 2024, 2025) both demonstrate that two massive accretion events are sufficient to generate a highly anisotropic stellar halo, with position- dependent break radii, power-law slopes, and large-scale north–south asymmetry. We refer the reader to J. A. S. Amarante et al. (2024), which presents a more detailed discussion on the position dependence of power indices. 6. Conclusions Based on the DR2 of the DESI MW Survey, we select a halo K giant sample in the range of 8 kpc < rGC < 200 kpc and perform detailed studies on the shape, orientation, radial 101 102 rq[kpc] 2 3 4 5 6 7 8 K giant BHB RRL MSTO Watkins 2009 Deason 2011 Xue 2015 Fukushima 2018 Medina 2018 Thomas 2018 Starkenburg 2019 Pieres 2020 Stringer 2021 Yang 2022 Han 2022 Yu 2024 Amarante 2024 Cavieres 2025 This work Figure 14. Radial density profiles of the MW stellar halo, comparing the result between this work (red solid line) and the literature. The vertical lines represent the break radii, and the horizontal lines denote the power indices within a given distance range, with the shaded range showing the uncertainties of power indices. The different line styles stand for various tracers. 16 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. density profile, and spatial anisotropies in the density distribution of the MW stellar halo, including local over- densities and the LMC density wake of our MW stellar halo. First, a forward triaxial ellipsoidal model with triple power- law functional form is fit to the observed spatial distribution of halo K giants after convolving with the survey selection functions. The best-fit model reveals two break radii at 16 and 76 kpc. We interpret the smaller break radius to be associated with the first pericentric passage of GSE and the larger break radius to be correlated with the infall of the LMC. The stellar halo major axis is tilted by 44° off the Galactic plane and 27° away from the Sun–Galactic center axis. With the same model, we find that the axis ratios are 1:p:q = 10:8:7, with the intermediate- and minor-axis lengths closer to each other than the major-axis length. By applying the forward model to halo K giants with different metallicities, we find that more metal-poor halo K giants have more extended radial profiles. We then allow the axis ratios and the major-axis orientation to vary with galactocentric distances, rGC. Interestingly, the inner stellar halo at rGC < ∼30 kpc is more oblate and more aligned with the disk, whereas the more distant stellar halo becomes more prolate and the major axis is more perpendicular to the MW disk. This indicates that the MW disk and inner stellar halo are misaligned with the outer stellar halo, while the outer stellar halo is more aligned with the VPOS of 11 classical MW satellite galaxies. Numerical simulations have predicted the change in the direction of angular momentum of galaxy disks and the inner halo, due to the infall of one or two massive satellites (e.g., F. A. Gómez et al. 2017), which may explain the twisted orientation and shape of the MW stellar halo in our measurement. Our MW stellar halo is full of substructures. We have successfully identified the HAC-N, HAC-S, and VOD over- densities associated with GSE. We identify a break radius at about 15 kpc in VOD and about 30 kpc in the HAC-N/S regions, consistent with J. J. Han et al. (2022). The perturbation of the LMC to the MW stellar halo is also investigated. The transient LMC density wake at larger distances due to dynamical friction is identified in the Pisces region, with a second break radius at ∼60 kpc. For fields free of local overdensities, we identify that the footprint in the northern Galactic cap is about 4.5 times more overdense than the footprint in the southern cap between 50 and 100 kpc. And there is a break radius at ∼100 kpc and an overdense peak at 90 kpc in the northern Galactic cap, which corresponds to a detection of the collective density wake of the LMC. The DESI data have allowed us to explore the stellar halo in greater depth and detail than ever before. Our findings show that the stellar halo is littered with substructures, and we are able to model and relate these overdensities to known accretion events and to explain the varying shape and density of halo stars across the sky. With the full 5 yr DESI sample, we anticipate a factor-of-2 reduction in the uncertainties of halo shape and orientation parameters, enabling a more precise measurement of the stellar halo twist signature and a better understanding of its origin. The measurements presented in this paper can be accessed in Zenodo at DOI: 10.5281/zenodo.18373975, which contains all data points for the figures presented in this work. Acknowledgments This work is supported by NSFC (12573022, 12595312, 12273021), the National Key R&D Program of China (2023YFA1605600, 2023YFA1605601), 111 project (No. B20019), and the Office of Science and Technology, Shanghai Municipal Government (grant Nos. 24DX1400100 and ZJ2023-ZD-001). We acknowledge the sponsorship from Yangyang Development Fund. The computations of this work are carried out on the National Energy Research Scientific Computing Center (NERSC) and the Gravity supercomputer at the Department of Astronomy, Shanghai Jiao Tong University. S.K. acknowledges support from the Science & Technology Facilities Council (STFC) grant ST/Y001001/1. M.V. acknowledges support from a NASA ATP award (80NSSC20K0509). A.P.C. acknowledges support from the Taiwan Ministry of Education Yushan Fellowship, MOE-113- YSFMS-0002-001-P2, and Taiwanese National Science and Technology Council grant 112-2112-M-007-009. L.B.e.S. acknowledges support from CNPq (Brazil) through a research productivity fellowship, grant No. [304873/2025-0]. We are grateful to Martin White for being our DESI publication handler. W.W. is grateful for useful discussions with Yanjun Sheng, Jianhui Lian, Baitian Tang, Jie Wang, Zhaozhou Li, Ling Zhu, Enci Wang, Hongxin Zhang, Yu Rong, Huiyuan Wang, and Bojun Tao. This material is based on work supported by the US Department of Energy (DOE), Office of Science, Office of High-Energy Physics, under contract No. DE– AC02–05CH11231, and by the National Energy Research Scientific Computing Center, a DOE Office of Science User Facility under the same contract. Additional support for DESI was provided by the US National Science Foundation (NSF), Division of Astronomical Sciences under contract No. AST- 0950945 to the NSF’s National Optical-Infrared Astronomy Research Laboratory; the Science and Technology Facilities Council of the United Kingdom; the Gordon and Betty Moore Foundation; the Heising-Simons Foundation; the French Alternative Energies and Atomic Energy Commission (CEA); the National Council of Humanities, Science and Technology of Mexico (CONAHCYT); the Ministry of Science, Innovation and Universities of Spain (MICIU/AEI/ 10.13039/501100011033); and the DESI Member Institutions: https://www.desi.lbl.gov/collaborating-institutions. Any opi- nions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the US National Science Foundation, the US Department of Energy, or any of the listed funding agencies. The authors are honored to be permitted to conduct scientific research on I’oligam Du’ag (Kitt Peak), a mountain with particular significance to the Tohono O’odham Nation. For the purpose of open access, the author has applied a Creative Commons Attribution (CC BY) licence to any Author Accepted Manuscript version arising from this submission. We thank the anonymous referee for his/her time and effort spent on reading and commenting on our paper. Appendix Effect of Distance Uncertainties In this appendix, we briefly discuss how the distance uncertainties of K giants could change the radial density 17 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. https://doi.org/10.5281/zenodo.18373975 https://www.desi.lbl.gov/collaborating-institutions profile and hence impact our results. We generate a mock density distribution with a triple power-law radial density profile without any observational uncertainties, which covers the whole sky (hereafter we name it the ideal profile). We assume that the parameters of the ideal profile are our best-fit values, as shown in Figure 3. We then perturb the ideal density distribution with the same distance error from the observation, and then we further convolve the perturbed density distribution with the DESI footprint and the angular + radial selection functions. We run MCMC again to get the best-fit parameters of the perturbed density distribution, as for real data. Figure 15 shows two radial density profiles of the ideal profile and the best-fit profile of the perturbed density distribution, with the green shaded region representing the model uncertainties of the perturbed density profile. The perturbed density profile shows larger break radii and different power indices. However, the deviation between the ideal profile and the perturbed profile is still within the model uncertainties. We have repeated this process multiple times to get different realizations of the perturbed profile. For all realizations, the difference is smaller than the statistical error. Hence, we believe that the distance uncertainty is not prominent in this work. ORCID iDs Songting Liaa https://orcid.org/0000-0002-6469-8263 Wenting Wangaa https://orcid.org/0000-0002-5762-7571 Sergey E. Koposovaa https://orcid.org/0000-0003-2644-135X João A. S. Amaranteaa https://orcid.org/0000-0002-7662-5475 Alis J. Deasonaa https://orcid.org/0000-0001-6146-2645 Nathan R. Sandfordaa https://orcid.org/0000-0002-7393-3595 Ting S. Liaa https://orcid.org/0000-0002-9110-6163 Gustavo E. Medinaaa https://orcid.org/0000-0003-0105-9576 Jiaxin Hanaa https://orcid.org/0000-0002-8010-6715 Monica Valluriaa https://orcid.org/0000-0002-6257-2341 Oleg Y. Gnedinaa https://orcid.org/0000-0001-9852-9954 Andrew P. Cooperaa https://orcid.org/0000-0001-8274-158X Leandro Beraldo e Silvaaa https://orcid.org/0000-0002- 0740-1507 Carlos Frenkaa https://orcid.org/0000-0002-2338-716X Raymond G. Carlbergaa https://orcid.org/0000-0002- 7667-0081 Mika Lambertaa https://orcid.org/0000-0002-2527-8899 Tian Qiuaa https://orcid.org/0000-0002-4900-2088 Jessica Nicole Aguilaraa https://orcid.org/0000-0003-0822-452X Steven Ahlenaa https://orcid.org/0000-0001-6098-7247 Davide Bianchiaa https://orcid.org/0000-0001-9712-0006 David Brooksaa https://orcid.org/0000-0002-8458-5047 Axel de la Macorraaa https://orcid.org/0000-0002-1769-1640 Peter Doelaa https://orcid.org/0000-0002-6397-4457 Jaime E. Forero-Romeroaa https://orcid.org/0000-0002- 2890-3725 Enrique Gaztañagaaa https://orcid.org/0000-0001-9632-0815 Satya Gontcho A Gontchoaa https://orcid.org/0000-0003- 3142-233X Gaston Gutierrezaa https://orcid.org/0000-0003-0825-0517 Robert Kehoeaa https://orcid.org/0000-0002-7101-697X Anthony Kreminaa https://orcid.org/0000-0001-6356-7424 Claire Lammanaa https://orcid.org/0000-0002-6731-9329 Martin Landriauaa https://orcid.org/0000-0003-1838-8528 Laurent Le Guillouaa https://orcid.org/0000-0001-7178-8868 Ramon Miquelaa https://orcid.org/0000-0002-6610-4836 Will Percivalaa https://orcid.org/0000-0002-0644-5727 Francisco Pradaaa https://orcid.org/0000-0001-7145-8674 Ignasi Pérez-Ràfolsaa https://orcid.org/0000-0001-6979-0125 Eusebio Sanchezaa https://orcid.org/0000-0002-9646-8198 David Schlegelaa https://orcid.org/0000-0002-5042-5088 Ray Sharplesaa https://orcid.org/0000-0003-3449-8583 Joseph Harry Silberaa https://orcid.org/0000-0002-3461-0320 David Sprayberryaa https://orcid.org/0000-0001-7583-6441 Gregory Tarléaa https://orcid.org/0000-0003-1704-0781 Hu Zouaa https://orcid.org/0000-0002-6684-3997 References Akaike, H. 1974, ITAC, 19, 716 Amarante, J. A. S., Debattista, V. P., Beraldo e Silva, L., Laporte, C. F. P., & Deg, N. 2022, ApJ, 937, 12 Amarante, J. A. S., Koposov, S. E., & Laporte, C. F. P. 2024, A&A, 690, A166 Baumgardt, H., & Vasiliev, E. 2021, MNRAS, 505, 5957 Bell, E. F., Zucker, D. B., Belokurov, V., et al. 2008, ApJ, 680, 295 Belokurov, V. 2013, NewAR, 57, 100 Belokurov, V., Deason, A. J., Erkal, D., et al. 2019, MNRAS, 488, L47 Belokurov, V., Erkal, D., Evans, N. W., Koposov, S. E., & Deason, A. J. 2018, MNRAS, 478, 611 Belokurov, V., Evans, N. W., Bell, E. F., et al. 2007, ApJL, 657, L89 Belokurov, V., Koposov, S. E., Evans, N. W., et al. 2013, MNRAS, 437, 116 Bett, P., Eke, V., Frenk, C. S., Jenkins, A., & Okamoto, T. 2010, MNRAS, 404, 1137 Bett, P. E., & Frenk, C. S. 2012, MNRAS, 420, 3324 Bett, P. E., & Frenk, C. S. 2016, MNRAS, 461, 1338 Bland-Hawthorn, J., & Gerhard, O. 2016, ARA&A, 54, 529 Bullock, J. S., & Johnston, K. V. 2005, ApJ, 635, 931 Byström, A., Koposov, S. E., Lilleengen, S., et al. 2025, MNRAS, 542, 560 Cautun, M., Bose, S., Frenk, C. S., et al. 2015a, MNRAS, 452, 3838 Cautun, M., Wang, W., Frenk, C. S., & Sawala, T. 2015b, MNRAS, 449, 2576 Cavieres, M., Chanamé, J., Navarrete, C., et al. 2025, ApJ, 983, 83 Chandra, V., Naidu, R. P., Conroy, C., et al. 2023a, ApJ, 951, 26 Chandra, V., Naidu, R. P., Conroy, C., et al. 2023b, ApJ, 956, 110 Chandra, V., Naidu, R. P., Conroy, C., et al. 2025, ApJ, 988, 156 Chua, K. T. E., Vogelsberger, M., Pillepich, A., & Hernquist, L. 2022, MNRAS, 515, 2681 Conroy, C., Bonaca, A., Cargile, P., et al. 2019, ApJ, 883, 107 Conroy, C., Naidu, R. P., Garavito-Camargo, N., et al. 2021, Natur, 592, 534 Cooper, A. P., Koposov, S. E., Allende Prieto, C., et al. 2023, ApJ, 947, 37 Cui, X.-Q., Zhao, Y.-H., Chu, Y.-Q., et al. 2012, RAA, 12, 1197 Deason, A. J., Belokurov, V., & Evans, N. W. 2011, MNRAS, 416, 2903 101 102 rq[kpc] 10 7 10 6 10 5 10 4 10 3 10 2 [k pc 3 ] Ideal Profile Perturbed Profile Figure 15. The black solid line refers to the selection-effect-free density profile based on a mock realization (ideal profile) of stars. The green line is the selection-effect-free best-fit model for the ideal profile perturbed with observational error. The green shaded region shows the statistical uncertainties of the best-fit model. The vertical dashed lines stand for the break radii of the original profile (black) and the recovered profile (green) after perturbation. Both profiles are renormalized to unity. 18 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. https://orcid.org/0000-0002-6469-8263 https://orcid.org/0000-0002-5762-7571 https://orcid.org/0000-0003-2644-135X https://orcid.org/0000-0002-7662-5475 https://orcid.org/0000-0001-6146-2645 https://orcid.org/0000-0002-7393-3595 https://orcid.org/0000-0002-9110-6163 https://orcid.org/0000-0003-0105-9576 https://orcid.org/0000-0002-8010-6715 https://orcid.org/0000-0002-6257-2341 https://orcid.org/0000-0001-9852-9954 https://orcid.org/0000-0001-8274-158X https://orcid.org/0000-0002-0740-1507 https://orcid.org/0000-0002-0740-1507 https://orcid.org/0000-0002-2338-716X https://orcid.org/0000-0002-7667-0081 https://orcid.org/0000-0002-7667-0081 https://orcid.org/0000-0002-2527-8899 https://orcid.org/0000-0002-4900-2088 https://orcid.org/0000-0003-0822-452X https://orcid.org/0000-0001-6098-7247 https://orcid.org/0000-0001-9712-0006 https://orcid.org/0000-0002-8458-5047 https://orcid.org/0000-0002-1769-1640 https://orcid.org/0000-0002-6397-4457 https://orcid.org/0000-0002-2890-3725 https://orcid.org/0000-0002-2890-3725 https://orcid.org/0000-0001-9632-0815 https://orcid.org/0000-0003-3142-233X https://orcid.org/0000-0003-3142-233X https://orcid.org/0000-0003-0825-0517 https://orcid.org/0000-0002-7101-697X https://orcid.org/0000-0001-6356-7424 https://orcid.org/0000-0002-6731-9329 https://orcid.org/0000-0003-1838-8528 https://orcid.org/0000-0001-7178-8868 https://orcid.org/0000-0002-6610-4836 https://orcid.org/0000-0002-0644-5727 https://orcid.org/0000-0001-7145-8674 https://orcid.org/0000-0001-6979-0125 https://orcid.org/0000-0002-9646-8198 https://orcid.org/0000-0002-5042-5088 https://orcid.org/0000-0003-3449-8583 https://orcid.org/0000-0002-3461-0320 https://orcid.org/0000-0001-7583-6441 https://orcid.org/0000-0003-1704-0781 https://orcid.org/0000-0002-6684-3997 https://doi.org/10.1109/TAC.1974.1100705 https://ui.adsabs.harvard.edu/abs/1974ITAC...19..716A/abstract https://doi.org/10.3847/1538-4357/ac8b0d https://ui.adsabs.harvard.edu/abs/2022ApJ...937...12A/abstract https://doi.org/10.1051/0004-6361/202450351 https://ui.adsabs.harvard.edu/abs/2024A&A...690A.166A/abstract https://doi.org/10.1093/mnras/stab1474 https://ui.adsabs.harvard.edu/abs/2021MNRAS.505.5957B/abstract https://doi.org/10.1086/588032 https://ui.adsabs.harvard.edu/abs/2008ApJ...680..295B/abstract https://doi.org/10.1016/j.newar.2013.07.001 https://ui.adsabs.harvard.edu/abs/2013NewAR..57..100B/abstract https://doi.org/10.1093/mnrasl/slz101 https://ui.adsabs.harvard.edu/abs/2019MNRAS.488L..47B/abstract https://doi.org/10.1093/mnras/sty982 https://ui.adsabs.harvard.edu/abs/2018MNRAS.478..611B/abstract https://doi.org/10.1086/513144 https://ui.adsabs.harvard.edu/abs/2007ApJ...657L..89B/abstract https://doi.org/10.1093/mnras/stt1862 https://ui.adsabs.harvard.edu/abs/2014MNRAS.437..116B/abstract https://doi.org/10.1111/j.1365-2966.2010.16368.x https://ui.adsabs.harvard.edu/abs/2010MNRAS.404.1137B/abstract https://ui.adsabs.harvard.edu/abs/2010MNRAS.404.1137B/abstract https://doi.org/10.1111/j.1365-2966.2011.20275.x https://ui.adsabs.harvard.edu/abs/2012MNRAS.420.3324B/abstract https://doi.org/10.1093/mnras/stw1395 https://ui.adsabs.harvard.edu/abs/2016MNRAS.461.1338B/abstract https://doi.org/10.1146/annurev-astro-081915-023441 https://ui.adsabs.harvard.edu/abs/2016ARA&A..54..529B/abstract https://doi.org/10.1086/497422 https://ui.adsabs.harvard.edu/abs/2005ApJ...635..931B/abstract https://doi.org/10.1093/mnras/staf1219 https://ui.adsabs.harvard.edu/abs/2025MNRAS.542..560B/abstract https://doi.org/10.1093/mnras/stv1557 https://ui.adsabs.harvard.edu/abs/2015MNRAS.452.3838C/abstract https://doi.org/10.1093/mnras/stv490 https://ui.adsabs.harvard.edu/abs/2015MNRAS.449.2576C/abstract https://doi.org/10.3847/1538-4357/adbf08 https://ui.adsabs.harvard.edu/abs/2025ApJ...983...83C/abstract https://doi.org/10.3847/1538-4357/accf13 https://ui.adsabs.harvard.edu/abs/2023ApJ...951...26C/abstract https://doi.org/10.3847/1538-4357/acf7bf https://ui.adsabs.harvard.edu/abs/2023ApJ...956..110C/abstract https://doi.org/10.3847/1538-4357/addab6 https://ui.adsabs.harvard.edu/abs/2025ApJ...988..156C/abstract https://doi.org/10.1093/mnras/stac1897 https://ui.adsabs.harvard.edu/abs/2022MNRAS.515.2681C/abstract https://doi.org/10.3847/1538-4357/ab38b8 https://ui.adsabs.harvard.edu/abs/2019ApJ...883..107C/abstract https://doi.org/10.1038/s41586-021-03385-7 https://ui.adsabs.harvard.edu/abs/2021Natur.592..534C/abstract https://doi.org/10.3847/1538-4357/acb3c0 https://ui.adsabs.harvard.edu/abs/2023ApJ...947...37C/abstract https://doi.org/10.1088/1674-4527/12/9/003 https://ui.adsabs.harvard.edu/abs/2012RAA....12.1197C/abstract https://doi.org/10.1111/j.1365-2966.2011.19237.x https://ui.adsabs.harvard.edu/abs/2011MNRAS.416.2903D/abstract Deason, A. J., Belokurov, V., Evans, N. W., & Johnston, K. V. 2013, ApJ, 763, 113 Deason, A. J., Belokurov, V., & Sanders, J. L. 2019, MNRAS, 490, 3426 DESI Collaboration, Abdul Karim, M., Aguilar, J., et al. 2025, PhRvD, 112, 083515 DESI Collaboration, Abareshi, B., Aguilar, J., et al. 2022, ApJ, 164, 207 DESI Collaboration, Adame, A. G., Aguilar, J., et al. 2024, AJ, 168, 58 DESI Collaboration, Aghamousa, A., Aguilar, J., et al. 2016a, arXiv:1611. 00036 DESI Collaboration, Aghamousa, A., Aguilar, J., et al. 2016b, arXiv:1611. 00037 DESI Collaboration, Adame, A. G., Aguilar, J., et al. 2025, JCAP, 2025, 028 Eggen, O. J., Lynden-Bell, D., & Sandage, A. R. 1962, ApJ, 136, 748 Emami, R., Genel, S., Hernquist, L., et al. 2021, ApJ, 913, 36 Faccioli, L., Smith, M. C., Yuan, H. B., et al. 2014, ApJ, 788, 105 Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 Freeman, K., & Bland-Hawthorn, J. 2002, ARA&A, 40, 487 Fukushima, T., Chiba, M., Homma, D., et al. 2018, PASJ, 70, 69 Fukushima, T., Chiba, M., Tanaka, M., et al. 2019, PASJ, 71, 72 Gaia Collaboration, Prusti, T., de Bruijne, J. H. J., et al. 2016, A&A, 595, A1 Garavito-Camargo, N., Besla, G., Laporte, C. F. P., et al. 2019, ApJ, 884, 51 Gnedin, O. Y., Brown, W. R., Geller, M. J., & Kenyon, S. J. 2010, ApJL, 720, L108 Gómez, F. A., Grand, R. J. J., Monachesi, A., et al. 2017, MNRAS, 472, 3722 Guy, J., Bailey, S., Kremin, A., et al. 2023, AJ, 165, 144 Hahn, C., Wilson, M. J., Ruiz-Macias, O., et al. 2023, AJ, 165, 253 Han, J., Wang, W., Cole, S., & Frenk, C. S. 2016a, MNRAS, 456, 1003 Han, J., Wang, W., Cole, S., & Frenk, C. S. 2016b, MNRAS, 456, 1017 Han, J. J., Conroy, C., & Hernquist, L. 2023, NatAs, 7, 1481 Han, J. J., Conroy, C., Johnson, B. D., et al. 2022, AJ, 164, 249 Hawkins, M. R. S. 1984, MNRAS, 206, 433 Hayden, M. R., Bovy, J., Holtzman, J. A., et al. 2015, ApJ, 808, 132 Haywood, M., Di Matteo, P., Lehnert, M. D., et al. 2018, ApJ, 863, 113 He, J., Wang, W., Li, Z., et al. 2024, ApJ, 976, 187 Helmi, A. 2020, ARA&A, 58, 205 Helmi, A., Babusiaux, C., Koppelman, H. H., et al. 2018, Natur, 563, 85 Helmi, A., White, S. D. M., de Zeeuw, P. T., & Zhao, H. 1999, Natur, 402, 53 Hernitschek, N., Cohen, J. G., Rix, H.-W., et al. 2018, ApJ, 859, 31 Huang, Y., Beers, T. C., Yuan, H.-B., et al. 2023, ApJ, 957, 65 Iorio, G., & Belokurov, V. 2019, MNRAS, 482, 3868 Johnston, K. V. 2016, ASSL, 420, 141 Jurić, M., Ivezić, Ž., Brooks, A., et al. 2008, ApJ, 673, 864 Kirk, D., Brown, M. L., Hoekstra, H., et al. 2015, SSRv, 193, 139 Kollmeier, J. A., Gould, A., Shectman, S., et al. 2009, ApJL, 705, L158 Koposov, S. E., Allende Prieto, C., Cooper, A. P., et al. 2024, MNRAS, 533, 1012 Koposov, S. E., Li, T. S., Prieto, C. A., et al. 2025, arXiv:2505.14787 Koppelman, H. H., Bos, R. O. Y., & Helmi, A. 2020, A&A, 642, L18 Lancaster, L., Koposov, S. E., Belokurov, V., Evans, N. W., & Deason, A. J. 2019, MNRAS, 486, 378 Lane, J. M. M., Bovy, J., & Mackereth, J. T. 2023, MNRAS, 526, 1209 Laporte, C. F. P., White, S. D. M., Naab, T., & Gao, L. 2013, MNRAS, 435, 901 Li, S., Wang, W., Koposov, S. E., et al. 2025, AJ, 170, 171 Li, Z., Han, J., Wang, W., et al. 2025, MNRAS, 538, 1442 Liu, C., Xu, Y., Wan, J.-C., et al. 2017, RAA, 17, 096 López-Corredoira, M., Tang, X. C., Tian, H., et al. 2024, A&A, 684, A135 Lowing, B., Wang, W., Cooper, A., et al. 2015, MNRAS, 446, 2274 Luo, A. L., Zhao, Y.-H., Zhao, G., et al. 2015, RAA, 15, 1095 Majewski, S. R., Schiavon, R. P., Frinchaboy, P. M., et al. 2017, AJ, 154, 94 Martell, S. L., Sharma, S., Buder, S., et al. 2017, MNRAS, 465, 3203 Medina, G. E., Li, T. S., Koposov, S. E., et al. 2025, arXiv:2504.02924 Medina, G. E., Muñoz, R. R., Carlin, J. L., et al. 2024, MNRAS, 531, 4762 Medina, G. E., Muñoz, R. R., Vivas, A. K., et al. 2018, ApJ, 855, 43 Miller, T. N., Doel, P., Gutierrez, G., et al. 2024, AJ, 168, 95 Morrison, H. L., Mateo, M., Olszewski, E. W., et al. 2000, AJ, 119, 2254 Naidu, R. P., Conroy, C., Bonaca, A., et al. 2021, ApJ, 923, 92 Newberg, H. J., Yanny, B., Rockosi, C., et al. 2002, ApJ, 569, 245 Nibauer, J., & Bonaca, A. 2025, ApJL, 985, L22 Nie, J. D., Smith, M. C., Belokurov, V., et al. 2015, ApJ, 810, 153 Okumura, T., & Jing, Y. P. 2009, ApJL, 694, L83 Okumura, T., Jing, Y. P., & Li, C. 2009, ApJ, 694, 214 Pace, A. B., Erkal, D., & Li, T. S. 2022, ApJ, 940, 136 Palau, C. G., & Miralda-Escudé, J. 2023, MNRAS, 524, 2124 Palau, C. G., Wang, W., & Han, J. 2025, MNRAS, 539, 2718 Pawlowski, M. S., & Kroupa, P. 2014, ApJ, 790, 74 Pawlowski, M. S., Pflamm-Altenburg, J., & Kroupa, P. 2012, MNRAS, 423, 1109 Perottoni, H. D., Limberg, G., Amarante, J. A. S., et al. 2022, ApJL, 936, L2 Pieres, A., Girardi, L., Balbinot, E., et al. 2020, MNRAS, 497, 1547 Pila-Díez, B., de Jong, J. T. A., Kuijken, K., van der Burg, R. F. J., & Hoekstra, H. 2015, A&A, 579, A38 Poppett, C., Tyas, L., Aguilar, J., et al. 2024, AJ, 168, 245 Preston, G. W., Shectman, S. A., & Beers, T. C. 1991, ApJ, 375, 121 Sawala, T., Cautun, M., Frenk, C., et al. 2023, NatAs, 7, 481 Schlafly, E. F., Kirkby, D., Schlegel, D. J., et al. 2023, AJ, 166, 259 Searle, L., & Zinn, R. 1978, ApJ, 225, 357 Sesar, B., Ivezić, Ž., Grammer, S. H., et al. 2010, ApJ, 708, 717 Sesar, B., Ivezić, Ž., Lupton, R. H., et al. 2007, AJ, 134, 2236 Sesar, B., Ivezić, Ž., Stuart, J. S., et al. 2013, AJ, 146, 21 Sesar, B., Jurić, M., & Ivezić, Ž. 2011, ApJ, 731, 4 Shao, S., Cautun, M., Deason, A., & Frenk, C. S. 2021, MNRAS, 504, 6033 Shao, S., Cautun, M., & Frenk, C. S. 2019, MNRAS, 488, 1166 Shao, S., Cautun, M., Frenk, C. S., et al. 2018, MNRAS, 476, 1796 Sheng, Y., Ting, Y.-S., & Xue, X.-X. 2025, MNRAS, 544, 2434 Sheng, Y., Ting, Y.-S., Xue, X.-X., Chang, J., & Tian, H. 2024, MNRAS, 534, 2694 Simion, I. T., Belokurov, V., & Koposov, S. E. 2019, MNRAS, 482, 921 Simion, I. T., Belokurov, V., Irwin, M., & Koposov, S. E. 2014, MNRAS, 440, 161 Sommer-Larsen, J. 1987, MNRAS, 227, 21P Starkenburg, E., Youakim, K., Martin, N., et al. 2019, MNRAS, 490, 5757 Stringer, K. M., Drlica-Wagner, A., Macri, L., et al. 2021, ApJ, 911, 109 Thomas, G. F., McConnachie, A. W., Ibata, R. A., et al. 2018, MNRAS, 481, 5223 Viswanathan, A., Byström, A., Starkenburg, E., et al. 2024, arXiv:2408.17250 Vivas, A. K., & Zinn, R. 2006, AJ, 132, 714 Vivas, A. K., Zinn, R., Andrews, P., et al. 2001, ApJL, 554, L33 Watkins, L. L., Evans, N. W., Belokurov, V., et al. 2009, MNRAS, 398, 1757 Wetterer, C. J., & McGraw, J. T. 1996, AJ, 112, 1046 Wheeler, A. J., Abruzzo, M. W., Casey, A. R., & Ness, M. K. 2023, AJ, 165, 68 White, S. D. M., & Frenk, C. S. 1991, ApJ, 379, 52 Wu, W., Zhao, G., Xue, X.-X., Bird, S. A., & Yang, C. 2022, ApJ, 924, 23 Xu, K., Jing, Y. P., & Gao, H. 2023a, ApJ, 954, 2 Xu, K., Jing, Y. P., & Zhao, D. 2023b, ApJ, 957, 45 Xu, Y., Liu, C., Xue, X.-X., et al. 2018, MNRAS, 473, 1244 Xue, X.-X., Ma, Z., Rix, H.-W., et al. 2014, ApJ, 784, 170 Xue, X.-X., Rix, H.-W., Ma, Z., et al. 2015, ApJ, 809, 144 Yang, C., Zhu, L., Tahmasebzadeh, B., Xue, X.-X., & Liu, C. 2022, AJ, 164, 241 Yanny, B., Rockosi, C., Newberg, H. J., et al. 2009, AJ, 137, 4377 Ye, D., Du, C., Deng, M., et al. 2024, MNRAS, 532, 2584 Yu, F., Li, T. S., Speagle, J. S., et al. 2024, ApJ, 975, 81 Zhu, L., Cai, R., Kang, X., et al. 2026, A&A, 706, A193 Zhu, L., Xue, X.-X., Mao, S., Yang, C., & Zhang, L. 2025, A&A, 703, A43 Zonca, A., Singer, L., Lenz, D., et al. 2019, JOSS, 4, 1298 19 The Astrophysical Journal, 999:108 (19pp), 2026 March 1 Li et al. https://doi.org/10.1088/0004-637X/763/2/113 https://ui.adsabs.harvard.edu/abs/2013ApJ...763..113D/abstract https://ui.adsabs.harvard.edu/abs/2013ApJ...763..113D/abstract https://doi.org/10.1093/mnras/stz2793 https://ui.adsabs.harvard.edu/abs/2019MNRAS.490.3426D/abstract https://doi.org/10.1103/tr6y-kpc6 https://ui.adsabs.harvard.edu/abs/2025PhRvD.112h3515A/abstract https://ui.adsabs.harvard.edu/abs/2025PhRvD.112h3515A/abstract https://doi.org/10.3847/1538-3881/ac882b https://ui.adsabs.harvard.edu/abs/2022AJ....164..207D/abstract https://doi.org/10.3847/1538-3881/ad3217 https://ui.adsabs.harvard.edu/abs/2024AJ....168...58D/abstract http://arXiv.org/abs/1611.00036 http://arXiv.org/abs/1611.00036 http://arXiv.org/abs/1611.00037 http://arXiv.org/abs/1611.00037 https://doi.org/10.1088/1475-7516/2025/07/028 https://ui.adsabs.harvard.edu/abs/2025JCAP...07..028A/abstract https://doi.org/10.1086/147433 https://ui.adsabs.harvard.edu/abs/1962ApJ...136..748E/abstract https://doi.org/10.3847/1538-4357/abf147 https://ui.adsabs.harvard.edu/abs/2021ApJ...913...36E/abstract https://doi.org/10.1088/0004-637X/788/2/105 https://ui.adsabs.harvard.edu/abs/2014ApJ...788..105F/abstract https://doi.org/10.1086/670067 https://ui.adsabs.harvard.edu/abs/2013PASP..125..306F/abstract https://ui.adsabs.harvard.edu/abs/2013PASP..125..306F/abstract https://doi.org/10.1146/annurev.astro.40.060401.093840 https://ui.adsabs.harvard.edu/abs/2002ARA&A..40..487F/abstract https://doi.org/10.1093/pasj/psy060 https://ui.adsabs.harvard.edu/abs/2018PASJ...70...69F/abstract https://doi.org/10.1093/pasj/psz052 https://ui.adsabs.harvard.edu/abs/2019PASJ...71...72F/abstract https://doi.org/10.1051/0004-6361/201629272 https://ui.adsabs.harvard.edu/abs/2016A&A...595A...1G/abstract https://doi.org/10.3847/1538-4357/ab32eb https://ui.adsabs.harvard.edu/abs/2019ApJ...884...51G/abstract https://doi.org/10.1088/2041-8205/720/1/L108 https://ui.adsabs.harvard.edu/abs/2010ApJ...720L.108G/abstract https://ui.adsabs.harvard.edu/abs/2010ApJ...720L.108G/abstract https://doi.org/10.1093/mnras/stx2149 https://ui.adsabs.harvard.edu/abs/2017MNRAS.472.3722G/abstract https://doi.org/10.3847/1538-3881/acb212 https://ui.adsabs.harvard.edu/abs/2023AJ....165..144G/abstract https://doi.org/10.3847/1538-3881/accff8 https://ui.adsabs.harvard.edu/abs/2023AJ....165..253H/abstract https://doi.org/10.1093/mnras/stv2707 https://ui.adsabs.harvard.edu/abs/2016MNRAS.456.1003H/abstract https://doi.org/10.1093/mnras/stv2522 https://ui.adsabs.harvard.edu/abs/2016MNRAS.456.1017H/abstract https://doi.org/10.1038/s41550-023-02076-9 https://ui.adsabs.harvard.edu/abs/2023NatAs...7.1481H/abstract https://doi.org/10.3847/1538-3881/ac97e9 https://ui.adsabs.harvard.edu/abs/2022AJ....164..249H/abstract https://doi.org/10.1093/mnras/206.3.433 https://ui.adsabs.harvard.edu/abs/1984MNRAS.206..433H/abstract https://doi.org/10.1088/0004-637X/808/2/132 https://ui.adsabs.harvard.edu/abs/2015ApJ...808..132H/abstract https://doi.org/10.3847/1538-4357/aad235 https://ui.adsabs.harvard.edu/abs/2018ApJ...863..113H/abstract https://doi.org/10.3847/1538-4357/ad8882 https://ui.adsabs.harvard.edu/abs/2024ApJ...976..187H/abstract https://doi.org/10.1146/annurev-astro-032620-021917 https://ui.adsabs.harvard.edu/abs/2020ARA&A..58..205H/abstract https://doi.org/10.1038/s41586-018-0625-x https://ui.adsabs.harvard.edu/abs/2018Natur.563...85H/abstract https://doi.org/10.1038/46980 https://ui.adsabs.harvard.edu/abs/1999Natur.402...53H/abstract https://doi.org/10.3847/1538-4357/aabfbb https://ui.adsabs.harvard.edu/abs/2018ApJ...859...31H/abstract https://doi.org/10.3847/1538-4357/ace628 https://ui.adsabs.harvard.edu/abs/2023ApJ...957...65H/abstract https://doi.org/10.1093/mnras/sty2806 https://ui.adsabs.harvard.edu/abs/2019MNRAS.482.3868I/abstract https://doi.org/10.1007/978-3-319-19336-6_6 https://ui.adsabs.harvard.edu/abs/2016ASSL..420..141J/abstract https://doi.org/10.1086/523619 https://ui.adsabs.harvard.edu/abs/2008ApJ...673..864J/abstract https://doi.org/10.1007/s11214-015-0213-4 https://ui.adsabs.harvard.edu/abs/2015SSRv..193..139K/abstract https://doi.org/10.1088/0004-637X/705/2/L158 https://ui.adsabs.harvard.edu/abs/2009ApJ...705L.158K/abstract https://doi.org/10.1093/mnras/stae1842 https://ui.adsabs.harvard.edu/abs/2024MNRAS.533.1012K/abstract https://ui.adsabs.harvard.edu/abs/2024MNRAS.533.1012K/abstract http://arXiv.org/abs/2505.14787 https://doi.org/10.1051/0004-6361/202038652 https://ui.adsabs.harvard.edu/abs/2020A&A...642L..18K/abstract https://doi.org/10.1093/mnras/stz853 https://ui.adsabs.harvard.edu/abs/2019MNRAS.486..378L/abstract https://doi.org/10.1093/mnras/stad2834 https://ui.adsabs.harvard.edu/abs/2023MNRAS.526.1209L/abstract https://doi.org/10.1093/mnras/stt912 https://ui.adsabs.harvard.edu/abs/2013MNRAS.435..901L/abstract https://ui.adsabs.harvard.edu/abs/2013MNRAS.435..901L/abstract https://doi.org/10.3847/1538-3881/adf1a0 https://ui.adsabs.harvard.edu/abs/2025AJ....170..171L/abstract https://doi.org/10.1093/mnras/staf358 https://ui.adsabs.harvard.edu/abs/2025MNRAS.538.1442L/abstract https://doi.org/10.1088/1674-4527/17/9/96 https://ui.adsabs.harvard.edu/abs/2017RAA....17...96L/abstract https://doi.org/10.1051/0004-6361/202348781 https://ui.adsabs.harvard.edu/abs/2024A&A...684A.135L/abstract https://doi.org/10.1093/mnras/stu2257 https://ui.adsabs.harvard.edu/abs/2015MNRAS.446.2274L/abstract https://doi.org/10.1088/1674-4527/15/8/002 https://ui.adsabs.harvard.edu/abs/2015RAA....15.1095L/abstract https://doi.org/10.3847/1538-3881/aa784d https://ui.adsabs.harvard.edu/abs/2017AJ....154...94M/abstract https://doi.org/10.1093/mnras/stw2835 https://ui.adsabs.harvard.edu/abs/2017MNRAS.465.3203M/abstract https://arxiv.org/abs/2504.02924 https://doi.org/10.1093/mnras/stae1137 https://ui.adsabs.harvard.edu/abs/2024MNRAS.531.4762M/abstract https://doi.org/10.3847/1538-4357/aaad02 https://ui.adsabs.harvard.edu/abs/2018ApJ...855...43M/abstract https://doi.org/10.3847/1538-3881/ad45fe https://ui.adsabs.harvard.edu/abs/2024AJ....168...95M/abstract https://doi.org/10.1086/301357 https://ui.adsabs.harvard.edu/abs/2000AJ....119.2254M/abstract https://doi.org/10.3847/1538-4357/ac2d2d https://ui.adsabs.harvard.edu/abs/2021ApJ...923...92N/abstract https://doi.org/10.1086/338983 https://ui.adsabs.harvard.edu/abs/2002ApJ...569..245N/abstract https://doi.org/10.3847/2041-8213/add0a9 https://ui.adsabs.harvard.edu/abs/2025ApJ...985L..22N/abstract https://doi.org/10.1088/0004-637X/810/2/153 https://ui.adsabs.harvard.edu/abs/2015ApJ...810..153N/abstract https://doi.org/10.1088/0004-637X/694/1/L83 https://ui.adsabs.harvard.edu/abs/2009ApJ...694L..83O/abstract https://doi.org/10.1088/0004-637X/694/1/214 https://ui.adsabs.harvard.edu/abs/2009ApJ...694..214O/abstract https://doi.org/10.3847/1538-4357/ac997b https://ui.adsabs.harvard.edu/abs/2022ApJ...940..136P/abstract https://doi.org/10.1093/mnras/stad1930 https://ui.adsabs.harvard.edu/abs/2023MNRAS.524.2124P/abstract https://doi.org/10.1093/mnras/staf658 https://ui.adsabs.harvard.edu/abs/2025MNRAS.539.2718P/abstract https://doi.org/10.1088/0004-637X/790/1/74 https://ui.adsabs.harvard.edu/abs/2014ApJ...790...74P/abstract https://doi.org/10.1111/j.1365-2966.2012.20937.x https://ui.adsabs.harvard.edu/abs/2012MNRAS.423.1109P/abstract https://ui.adsabs.harvard.edu/abs/2012MNRAS.423.1109P/abstract https://doi.org/10.3847/2041-8213/ac88d6 https://ui.adsabs.harvard.edu/abs/2022ApJ...936L...2P/abstract https://doi.org/10.1093/mnras/staa1980 https://ui.adsabs.harvard.edu/abs/2020MNRAS.497.1547P/abstract https://doi.org/10.1051/0004-6361/201425457 https://ui.adsabs.harvard.edu/abs/2015A&A...579A..38P/abstract https://doi.org/10.3847/1538-3881/ad76a4 https://ui.adsabs.harvard.edu/abs/2024AJ....168..245P/abstract https://doi.org/10.1086/170175 https://ui.adsabs.harvard.edu/abs/1991ApJ...375..121P/abstract https://doi.org/10.1038/s41550-022-01856-z https://ui.adsabs.harvard.edu/abs/2023NatAs...7..481S/abstract https://doi.org/10.3847/1538-3881/ad0832 https://ui.adsabs.harvard.edu/abs/2023AJ....166..259S/abstract https://doi.org/10.1086/156499 https://ui.adsabs.harvard.edu/abs/1978ApJ...225..357S/abstract https://doi.org/10.1088/0004-637X/708/1/717 https://ui.adsabs.harvard.edu/abs/2010ApJ...708..717S/abstract https://doi.org/10.1086/521819 https://ui.adsabs.harvard.edu/abs/2007AJ....134.2236S/abstract https://doi.org/10.1088/0004-6256/146/2/21 https://ui.adsabs.harvard.edu/abs/2013AJ....146...21S/abstract https://doi.org/10.1088/0004-637X/731/1/4 https://ui.adsabs.harvard.edu/abs/2011ApJ...731....4S/abstract https://doi.org/10.1093/mnras/staa3883 https://ui.adsabs.harvard.edu/abs/2021MNRAS.504.6033S/abstract https://doi.org/10.1093/mnras/stz1741 https://ui.adsabs.harvard.edu/abs/2019MNRAS.488.1166S/abstract https://doi.org/10.1093/mnras/sty343 https://ui.adsabs.harvard.edu/abs/2018MNRAS.476.1796S/abstract https://doi.org/10.1093/mnras/staf1780 https://ui.adsabs.harvard.edu/abs/2025MNRAS.544.2434S/abstract https://doi.org/10.1093/mnras/stae2259 https://ui.adsabs.harvard.edu/abs/2024MNRAS.534.2694S/abstract https://ui.adsabs.harvard.edu/abs/2024MNRAS.534.2694S/abstract https://doi.org/10.1093/mnras/sty2744 https://ui.adsabs.harvard.edu/abs/2019MNRAS.482..921S/abstract https://doi.org/10.1093/mnras/stu133 https://ui.adsabs.harvard.edu/abs/2014MNRAS.440..161S/abstract https://ui.adsabs.harvard.edu/abs/2014MNRAS.440..161S/abstract https://doi.org/10.1093/mnras/227.1.21P https://ui.adsabs.harvard.edu/abs/1987MNRAS.227P..21S/abstract https://doi.org/10.1093/mnras/stz2935 https://ui.adsabs.harvard.edu/abs/2019MNRAS.490.5757S/abstract https://doi.org/10.3847/1538-4357/abe873 https://ui.adsabs.harvard.edu/abs/2021ApJ...911..109S/abstract https://doi.org/10.1093/mnras/sty2604 https://ui.adsabs.harvard.edu/abs/2018MNRAS.481.5223T/abstract https://ui.adsabs.harvard.edu/abs/2018MNRAS.481.5223T/abstract https://arxiv.org/abs/2408.17250 https://doi.org/10.1086/505200 https://ui.adsabs.harvard.edu/abs/2006AJ....132..714V/abstract https://doi.org/10.1086/320915 https://ui.adsabs.harvard.edu/abs/2001ApJ...554L..33V/abstract https://doi.org/10.1111/j.1365-2966.2009.15242.x https://ui.adsabs.harvard.edu/abs/2009MNRAS.398.1757W/abstract https://doi.org/10.1086/118076 https://ui.adsabs.harvard.edu/abs/1996AJ....112.1046W/abstract https://doi.org/10.3847/1538-3881/acaaad https://ui.adsabs.harvard.edu/abs/2023AJ....165...68W/abstract https://ui.adsabs.harvard.edu/abs/2023AJ....165...68W/abstract https://doi.org/10.1086/170483 https://ui.adsabs.harvard.edu/abs/1991ApJ...379...52W/abstract https://doi.org/10.3847/1538-4357/ac31ac https://ui.adsabs.harvard.edu/abs/2022ApJ...924...23W/abstract https://doi.org/10.3847/1538-4357/ace62b https://ui.adsabs.harvard.edu/abs/2023ApJ...954....2X/abstract https://doi.org/10.3847/1538-4357/acf835 https://ui.adsabs.harvard.edu/abs/2023ApJ...957...45X/abstract https://doi.org/10.1093/mnras/stx2361 https://ui.adsabs.harvard.edu/abs/2018MNRAS.473.1244X/abstract https://doi.org/10.1088/0004-637X/784/2/170 https://ui.adsabs.harvard.edu/abs/2014ApJ...784..170X/abstract https://doi.org/10.1088/0004-637X/809/2/144 https://ui.adsabs.harvard.edu/abs/2015ApJ...809..144X/abstract https://doi.org/10.3847/1538-3881/ac9900 https://ui.adsabs.harvard.edu/abs/2022AJ....164..241Y/abstract https://ui.adsabs.harvard.edu/abs/2022AJ....164..241Y/abstract https://doi.org/10.1088/0004-6256/137/5/4377 https://ui.adsabs.harvard.edu/abs/2009AJ....137.4377Y/abstract https://doi.org/10.1093/mnras/stae1655 https://ui.adsabs.harvard.edu/abs/2024MNRAS.532.2584Y/abstract https://doi.org/10.3847/1538-4357/ad738f https://ui.adsabs.harvard.edu/abs/2024ApJ...975...81Y/abstract https://doi.org/10.1051/0004-6361/202557613 https://ui.adsabs.harvard.edu/abs/2026A&A...706A.193Z/abstract https://doi.org/10.1051/0004-6361/202556036 https://ui.adsabs.harvard.edu/abs/2025A&A...703A..43Z/abstract https://doi.org/10.21105/joss.01298 https://ui.adsabs.harvard.edu/abs/2019JOSS....4.1298Z/abstract 1. Introduction 2. Data 2.1. The DESI Milky Way Survey 2.2. MWS Stellar Parameter Catalog 2.3. Selection of K Giants 3. Methods 3.1. Correction of Incompleteness 3.1.1. Correcting Incompleteness of the Spectroscopic Survey 3.1.2. Correcting Incompleteness of Photometric Targets 3.2. Constructing the Stellar Halo with Forward Modeling 3.3. Likelihood Function 4. Results 4.1. Best-fit Triple Power-law Model 4.2. Radially Variable Model with Fixed Radial Density Profile 4.3. Spatial Anisotropies of Stellar Halo Density Distribution and the LMC Density Wake 4.4. Model Dependence on Metallicity 5. Discussions 5.1. Physical Interpretations on the Twisted Inner and Outer Stellar Halos 5.2. Comparison with Other Studies 6. Conclusions Appendix. Effect of Distance Uncertainties References