The remote response of the South Asian Monsoon to reduced dust emissions and Sahara greening during the middle Holocene

. Previous studies based on multiple paleoclimate archives suggested a prominent intensification of the South Asian Monsoon (SAM) during the mid-Holocene (MH, ~ 6000 years before present day). The main forcing that contributed to this intensification is related to changes in the Earth’s 20 orbital parameters. However, other key factors likely played important roles, including remote changes in vegetation cover and airborne dust emission. In particular, northern Africa also experienced much wetter conditions and a more mesic landscape than today during the MH (the so-called African Humid Period), leading to a large decrease in airborne dust globally. However, most modelling studies investigating the SAM changes during the Holocene overlooked the potential impacts of the vegetation 25 and dust emission changes that took place over northern Africa. Here, we use a set of simulations for the MH climate, in which vegetation over the Sahara and reduced dust concentrations are considered. Our results show that SAM rainfall is strongly affected by Saharan vegetation and dust concentrations, with a large increase in particular over northwestern India and a lengthening of the monsoon season. We propose that this remote influence is mediated by anomalies in Indian Ocean sea-surface temperatures 30 and may have shaped the evolution of the SAM during the termination of the African Humid Period.


Introduction
The South Asian Monsoon (SAM) directly affects the climate of the Indian subcontinent and indirectly influences far-afield regions through atmospheric and oceanic teleconnections [e.g., Lau, 1992;Liu et al., 2004]. Due to its key role for regional and global hydrological cycles, much attention has been 35 devoted to better understand and predict its variability on multiple timescales, including its long-term future changes (e.g., Huo and Peltier, 2020;Swapna et al., 2018). However, SAM future projections are highly uncertain (e.g., Huang et al., 2020), and even representing the recent trend and identifying its drivers has been challenging (e.g., Mishra et al., 2018) due to the relatively short modern observational record that spans roughly a century. Hence, investigating past SAM changes is of utmost importance to 40 better understand its dynamics and future evolution.
Dramatic shifts in the intensity of the SAM occurred at the end of the deglaciation (Bird et al., 2014;Campo et al., 1982;Dallmeyer et al., 2013;Fleitmann et al., 2003;Gill et al., 2017;Saraswat et al., 2013) when stronger boreal summer insolation, higher greenhouse gas concentrations, and shrinking ice 45 sheets triggered a strengthening of the northern hemisphere summer monsoon systems (Jalihal et al., 2019a;Sun et al., 2019). In particular, the increased orbital forcing enhanced moisture transport from the Indian Ocean to the Indian subcontinent, leading to increased monsoonal precipitation there (e.g., Dallmeyer et al., 2013;Texier et al., 2000). These changes occurred in parallel with a prolonged period of intense precipitation over north-western Africa -labelled the African Humid Period. The African 50 Humid Period spanned the early and middle Holocene (15,000 -4,000 years BP), and had far-reaching local and global climatic influences (Muschitiello et al., 2015;Pausata et al., 2017aPausata et al., , 2017bPiao et al., 2020;Sun et al., 2019). Locally, it coincided with a major intensification of the West African Monsoon (WAM) and a greening of the present-day Sahara Desert. Amongst its many remote impacts, the WAM strengthening contributed to the greening of the arid and semi-arid regions of east and south Asia (see 55 for a recent review ). Indeed, the large circulation changes instigated by the African Humid Period greening of the Sahara, together with the associated changes in sea surface temperatures, have likely complemented orbital changes in modulating the SAM (see (Texier et al., 2000) relative to land-surface changes alone).
[2017] and Messori et al. [2019] have shown that atmospheric dust loading profoundly affects monsoonal dynamics and, while changes in vegetation do lead to increased monsoonal precipitation, a better agreement with proxy data is only reached when a dust reduction is also simulated (see also Tierney et al., 2017). The latter strengthens the effects of land-surface changes, leading to a further increase and northward extension of the WAM. Egerer et al. [2018] have also shown that accounting 70 for both vegetation and dust feedbacks leads to a better match between model simulations and paleoclimate reconstructions from the north-western African margin. Another recent study (Thompson et al., 2019) has suggested a contribution from dust aerosol reduction of about 15-20% to the total rainfall over the Sahara, although these numbers may be dependent on the modelled dust optical properties and particle size range (Hopcroft and Valdes, 2019). 75 Through a set of sensitivity experiments performed with an Earth System Model, Pausata et al. [2017aPausata et al. [ , 2017b have shown that the strengthening of the WAM and the associated vegetation and dust feedbacks during the MH are able to affect the El Niño Southern Oscillation variability as well as tropical storm activity worldwide. Using the same set of simulations, Piao et al. [2020] show that a 80 vegetated Sahara leads to an enhancement of the western Pacific subtropical high, which in turn strengthens the East Asian Summer Monsoon. Sun et al. (2019) highlighted that Northern Hemisphere land monsoon precipitation significantly increases by over 30% under the effect of the Green Sahara.
However, a systematic evaluation of the joint impacts of atmospheric dust loading reductions, Saharan land-cover changes, and insolation changes on the SAM during the middle Holocene is lacking in 85 current literature.
https://doi.org/10.5194/cp-2020-142 Preprint. Discussion started: 21 November 2020 c Author(s) 2020. CC BY 4.0 License. 4 Here, we address this gap with the aim of providing insights into future SAM changes. Indeed, a number of recent studies have projected future increases in Sahelian precipitation (Biasutti, 2013; Giannini and Kaplan, 2019) associated with a surface greening and reduced dust emissions (Evan et al.,90 2016).
The remainder of the manuscript is organised as follows: The climate model used and the experimental design are described in Section 2. Next, we examine SAM changes during the summer, both at the surface and aloft (Section 3). A discussion and conclusions follow in Section 4. 95

Model description and experimental set-up
The study is based on a set of simulations performed with the fully coupled global climate model EC-Earth version 3.1. EC-Earth version 2 participated in the fifth phase of the Coupled Model Intercomparison Project (CMIP5) and version 3 will participate in CMIP6 (http://www.ec-earth.org/).
The model is comprised of the Integrated Forecasting System (IFS cycle 36r4) for the atmosphere, the 100 Nucleus for European Modelling of the Ocean version 2 (NEMO2) -for the ocean, and the Louvain-la-Neuve sea-ice Model version 3 (LIM3) for the sea-ice (Hazeleger et al., 2010;Yepes-Arbós et al., 2016). The IFS model includes the H-TESSEL land surface scheme and is run at T159 horizontal spectral resolution corresponding to roughly 1.125° in longitude and latitude, with 62 vertical levels.
NEMO has a 1° horizontal resolution except at the equator, where it increases to 1/3° (Sterl et al.,105 2012), and 46 vertical levels. The different components are coupled via the OASIS3 coupler. Relevant for this study, vegetation cover and monthly aerosol concentrations are prescribed in the model; however, the indirect effect of aerosols on clouds is not considered.
We analyse an MH experiment (MH PMIP ), which follows the protocol for the standard mid-Holocene 110 simulations in accordance with the third phase of the Paleoclimate Modelling Intercomparison Project (PMIP3) (Taylor et al., 2009(Taylor et al., , 2012 and three sensitivity experiments performed by Pausata et al. (2016) and Gaetani et al. (2017) (Table 1). The MH PMIP includes mid-Holocene orbital forcing and greenhouse gas concentrations, pre-industrial land cover, and airborne dust concentrations. The three https://doi.org/10.5194/cp-2020-142 Preprint. Discussion started: 21 November 2020 c Author(s) 2020. CC BY 4.0 License. 5 sensitivity experiments were carried out to investigate the effects of changes to land-cover conditions 115 and dust concentration in isolation as well as in combination. In the MH GS ('Green Sahara') setup, the vegetation type (and related parameters, see below) over the Sahara (defined as the area 11° -33°N, 15°W -35°E) is prescribed to be evergreen shrub, representing an idealised African Humid Period scenario, while dust concentration is left unaltered at its pre-industrial (PI) amounts. In the MH RD ('Reduced Dust') setup, the dust concentration over North Africa is reduced by up to 80% relative to 120 pre-industrial values (see Figs. 1 and S1 in Gaetani et al. (2017)), while the land-surface properties are kept to PI values. The final experiment (MH GS+RD ) considers the case where both vegetation and dust changes described above are simultaneously prescribed. The imposed changes in vegetation type correspond to important changes in surface albedo and leaf area index (LAI) as summarised in Table 2.
Under the Green Sahara scenario, the albedo decreases from 0.3 to 0.15, while the LAI increases from 125 0.2 to 2.6. For a more detailed description of these simulations, the reader is referred to Pausata et al. (2016) and Gaetani et al. (2017). These sensitivity experiments are compared to a PI simulation to investigate the role of each forcing in altering the SAM. The analysis focuses on the June-September (JJAS) period using the last 50 years of each experiment. Finally, the statistical significance of the differences between experiments at the 5% level is evaluated by a two-tailed Student's t test.

Results
This section discusses the SAM response in terms of local (Section 3.1) and large-scale changes (Section 3.2) to each forcing independently and together: orbital (MH PMIP ), orbital forcing and Sahara greening (MH GS ), orbital forcing and dust reduction (MH RD ), and orbital forcing, Sahara greening, and dust reduction (MH GS+RD ). In Section 3.3 we then compare the model findings to paleoclimate archives.

Precipitation
In the PI experiment, the SAM displays the most intense summertime (particularly June and July) precipitation over the west coast of the Indian subcontinent and the Himalayan foothills (Fig. 1a), in overall agreement with observations (Figs. A1 and A2). The MH PMIP experiment simulates a general 6 increase in SAM rainfall over South Asia compared to PI (Fig. 2a) as also shown by other PMIP model experiments (e.g., (Zhao and Harrison, 2012)), particularly over southern India and the Himalayan foothills. In contrast, decreased precipitation is seen over most of the Bay of Bengal, South China Sea, Indochina, and Thailand. This decrease in precipitation is a result of the reduced surface latent heat flux over the ocean as shown in (Jalihal et al., 2019b). This results in a decrease in the net energy flux into 145 the atmosphere over these regions, leading to a decline in precipitation. A precipitation anomaly dipole is simulated along the equatorial Indian Ocean, with increased precipitation to the west and decrease to the east. The greening of the Sahara (MH GS ) leads to a general intensification of the anomaly pattern simulated when only including orbital forcing (MH PMIP ; Fig. A3a). However, some peculiar characteristics emerge: in particular, the precipitation increases over a broad swathe of north-western 150 India and Pakistan, while it decreases over the Western Ghats (cf. panels a and b in Figure 2 and see also Figure A3a). The positive rainfall anomaly over the western equatorial Indian Ocean extends eastward, strongly reducing the negative precipitation anomaly in the eastern side of the basin. The reduction in precipitation over the Bay of Bengal, Indochina, and Thailand further intensifies. The reduced Saharan dust (MH RD ) leads to a pattern that is very similar -albeit with weaker anomalies -to 155 the orbital only forcing (MH PMIP ; Fig. A3b); however, the precipitation increase over southern India is confined to east of the Western Ghats, while a small decrease in rainfall is simulated along the western coast of the Indian subcontinent (cf. panels a and c in Figure 2 and see also Figure A3b Ghats is further enhanced in the MH GS+RD , while the increase over the Himalayan foothills is reduced compared to the MH GS , which is due to the effect of the dust reduction (cf. panels a and c in Figure 2 and see also Figure A3b). 165 While seasonal-mean precipitation determines the overall amount of water supplied, sub-seasonal changes in the monsoon, such as a shift in the onset and/or the withdrawal, are key to determining the length and hence the precipitation rate over the monsoonal season. The SAM in the PI simulation starts https://doi.org/10.5194/cp-2020-142 Preprint. Discussion started: 21 November 2020 c Author(s) 2020. CC BY 4.0 License. 7 in late May (Figs. 3 and A2a), with the monsoon then developing until early August and retreating in early September (Figs. 3 and A2a). In the MH PMIP experiment, the model simulates a delayed onset 170 south of 15°N, but not at higher latitudes (Fig. 3a). The withdrawal is, however, delayed at all latitudes, lengthening the overall duration of the monsoon by about one month. The Sahara greening (MH GS ) leads to a further lengthening of the monsoon season from April to October, and increased cumulative precipitation over a large part of the SAM region (Fig. 3b). The delayed onset is confined to the region well south of 10°N. While showing a lengthening of the monsoon season, the regions around the 15°N 175 latitude band show a decrease in precipitation between June and mid-August (Fig. 3b). Dust reduction (MH RD ) leads to a much stronger delay of the monsoon onset that extends up to 25°N compared to the MH PMIP experiment (Fig. 3c). The withdrawal of the monsoon is also delayed and resembles the MH PMIP simulation (Fig. 3c). In the MH GS+RD case, the distribution of rainfall is dominated by the Sahara greening, but the footprint of the dust reduction is visible at the lower latitudes. These display a 180 stronger decrease in precipitation than the MH GS simulation, in particular during the core monsoonal season -June to late August (Fig. 3d). Therefore, reduced dust seems to primarily reduce the South Asian monsoonal precipitation at low latitudes, while the greening of the Sahara increases the precipitation further north. 185

Surface temperature
The highest PI surface temperatures on the Indian subcontinent are simulated over its northwestern part, in agreement with observations ( Fig. A1). However, while the temperature pattern is similar to observations, our simulation displays a cold anomaly over a large part of the domain (  Figure 4 and see also Figure A5a). However, the warming is more pronounced than in MH PMIP over the bulk of the domain, with north-western India being an exception (Fig. A5a).
This may be linked to the simulated increase in rainfall in the region (Fig. 2b, A3a). The cold SST

Evapotranspiration
From an impacts-based perspective, changes in precipitation are only one part of the hydrological cycle, which also includes evaporation and, over land, transpiration as well. Therefore, it is important to also 215 investigate the evapotranspiration changes during the MH climate to better understand the impacts of changes in orbital forcing and Saharan vegetation and dust on the water budget of South Asia. In the PI experiment, weak evapotranspiration is simulated over the dry sub-tropical desert regions, while rates in excess of 3 mm/day are present over the Indian subcontinent ( Fig. 1c). In the MH PMIP experiment, the evapotranspiration is increased over north-western and south-eastern India, the Tibetan Plateau, and 220 Indochina (Fig. 5a). These regions are characterized by increases in precipitation and/or in temperature (Figs. 2 and 4), which enhance the evapotranspiration. On the contrary, the Indian Ocean displays a widespread decrease in evaporation rates except along the coast of Somalia. The Sahara greening (MH GS ) leads to a widespread increase in evapotranspiration across most of the Indian subcontinent, https://doi.org/10.5194/cp-2020-142 Preprint. Discussion started: 21 November 2020 c Author(s) 2020. CC BY 4.0 License. 9 enhancing the anomaly pattern simulated in the MH PMIP experiment (panels a and b in Figure 5 and see 225 also Figure A6a). The reduction in airborne dust (MH RD ) does not notably alter the evaporation over land compared to the orbital forcing only experiment (MH PMIP ); however, it significantly increases the evaporation over the Arabian Sea (panels a and c in Figure 5 and see also Figure A6b), due to the increase in incoming solar radiation. Finally, the combined forcing (MH GS+RD ) leads to mainly positive anomalies over land, as in the MH GS case, while the effects of dust reduction dominate over the Arabian 230 Sea and western Indian Ocean (Fig. 5d).

Changes in the large-scale monsoonal circulation
The PI sea-level pressure (SLP) pattern displays a thermal low over the Arabian Peninsula extending into the northern part of the Indian subcontinent ( Fig. 1d). This is associated with an anticyclonic circulation over the Indian Ocean leading to a strong westerly flow across the Indian subcontinent and 235 Indochina, which brings large amounts of moisture to these regions (Figs. 1a and d, A8a). The strong westerlies over the Arabian Sea favour upwelling and explain the origin of the "cold pool" in that region increases evaporation in the MH PMIP compared to the PI experiment (Fig. 5a). Conversely, the decreased evaporation in the northern Arabian Sea may be ascribed to the strengthened monsoonal flow, which 250 increases upwelling and in turn cools the region (Fig. 4a). Finally, the weakened westerly flow around https://doi.org/10.5194/cp-2020-142 Preprint. Discussion started: 21 November 2020 c Author(s) 2020. CC BY 4.0 License.
the southern tip of India may be responsible for decreased evaporation and a consequent increase in SSTs of that region. Under the Green Sahara conditions (MH GS ), the SLP anomaly pattern intensifies relative to the MH PMIP and shifts to the northwest, thus weakening the south-westerlies over the Arabian Sea, while strengthening the easterlies over the southern tip of India (panels a and b in Figure 6 and see 255 also Figure A7a). The latter anomaly can explain the decrease in precipitation over the western slopes of the Western Ghats and the increase on their eastern side. Although the south-westerly flow over the Arabian Sea is less intense than in the MH PMIP (Fig. A7a), the moisture advection is enhanced (Figs. 6b and A9a), which explains the increased precipitation and evapotranspiration over most of India (Figs. 2b, 5b, A6a). Indeed, the weakened atmospheric flow decreases the upwelling and in turn increases 260 SSTs, favouring more evaporation over the Arabian Sea (Fig. A6a). Reduced Saharan dust (MH RD ) results in a northward expansion of the Mascarene High in the southern Indian Ocean and a weakening of the Saudi Arabian heat low relative to MH PMIP experiment (panels a and c in Figure 6 and see also Figure A7b). This leads to a weakening of the Somali Jet, a weaker coastal upwelling in the Arabian Sea favouring modest warm SST anomalies there (Fig. 4c), and ultimately a weaker moisture transport 265 from the Arabian Sea to the southern half of the Indian subcontinent (panels a and c in Figure 7 and see also Figure A9b). The weakened low-level winds relative to MH PMIP are consistent with the significant decrease in precipitation over western India (Fig. A3b). Further east, there is a strengthened northwesterly flow over the Bay of Bengal extending towards the equatorial western Pacific, associated with a decreased moisture convergence over Bangladesh and north-eastern India relative to the MH PMIP 270 simulation (Fig. A9b). This circulation change causes a precipitation increase in the MH RD that is smaller than in the MH PMIP relative to the PI (Figs. 2c and A3b). When combining both Sahara greening and dust reduction (MH GS+RD ), SLP anomalies are mostly a linear combination of the two forcings ( Fig.   6d). In particular, the cyclonic footprint over the Indian Ocean and the easterly moisture transport from the Pacific to the Indian Ocean are both features of the MH GS experiment (Figs. A7 and A9). On the 275 other hand, over the Arabian Sea, both forcings contribute to a weakened westerly flow, albeit at slightly different latitudes.
https://doi.org/10.5194/cp-2020-142 Preprint. Discussion started: 21 November 2020 c Author(s) 2020. CC BY 4.0 License. 11 We next analyse the mid and upper-level circulation associated with the monsoonal flows. The PI 500-hPa vertical velocity field shows a strong ascending flow across the tropics during the monsoon season 280 (Fig. A8b), matching the areas of low SLP shown in Figure 1d, with the clear exception of the areas under thermal low pressures (e.g., Saudi Arabia and Iran). Subsidence is largely limited to the west Arabian Sea and Somali peninsula (Fig. A8b). Additionally, strong subsidence occurs over the desert regions of the Arabian Peninsula and Iran. Changes in orbital forcing (MH PMIP ) drive a strengthened upward motion over the western north-equatorial Indian Ocean, southern India, and Himalayan foothills 285 (Fig. 8a). This favours cloud formation and is consistent with increased precipitation over these regions  (Fig. A10b). The result of the Sahara greening and dust reduction forcing (MH GS+RD ) over Asia is to a great extent a linear combination of the two separate forcings (Fig. 8d), as was indeed the case for the other variables analysed here.
We next discuss the upper-level velocity potential and divergent winds, which provide a framework to 300 analyse the regional anomalies in the context of the large-scale tropical overturning circulation. The PI experiment shows a divergent flow emanating from Southeast Asia towards the surrounding Asian Monsoon regions (contour lines in Fig. 9), which is consistent with the low SLP there (Fig. 1d). In the MH PMIP , the whole pattern of velocity potential and the centres of divergence/convergence are shifted westward (Fig. 9a), with dipole anomalies centred over Northern Africa and the Arabian Peninsula 305 (negative velocity potential/divergence) and South America (positive velocity potential/convergence).

12
The divergence over the north-western Indian subcontinent is strengthened, which implies an intensified low-level convergence and hence stronger precipitation in the region. The greening of the Sahara (MH GS ) further intensifies the anomaly pattern seen in the MH PMIP experiment (cf. panels a and b in Figure 9 and see also Figure A11a). The dust reduction experiment contributes to a strong positive 310 anomaly in velocity potential over the Arabian Sea relative to MH PMIP (cf. panels a and c in Figure 9 and see also Figure A11b), thus weakening upper tropospheric divergence and the lower tropospheric convergence. The Green Sahara-reduced dust (MH GS+RD ) experiment resembles the MH GS forcing, but the anomalies are reduced due to the effect of dust reduction ( Fig. 9 and A11c). 315 The anomalies in velocity potential are negative over both India and the Bay of Bengal, albeit with smaller magnitudes over the latter region. Therefore, the decrease in precipitation over the Bay of This further leads to a decrease in the net energy flux into the atmosphere (top + bottom) over the Bay 325 of Bengal. Since precipitation is proportional to the net energy flux into the atmosphere, precipitation over the Bay of Bengal decreases (Jalihal et al., 2019b). In general, as the precipitation over the WEIO and NEA increases, there is a corresponding decrease in latent heat flux over the Bay of Bengal. MH GS shows the largest increase in precipitation over the WEIO and NEA (Fig. 10a). Proportionately, the decrease in latent heat flux over the Bay of Bengal is also the largest. On the other hand, the weakest 330 increase in precipitation over the WEIO and NEA regions is simulated in the MH RD . The associated reduction in latent heat flux over the Bay of Bengal is also the smallest. As the latent heat flux decreases, it leads to a larger reduction in precipitation over the Bay of Bengal (Fig. 10b). This change in latent heat flux is due to the impact of precipitation over the WEIO and NEA on wind speed over the https://doi.org/10.5194/cp-2020-142 Preprint. Discussion started: 21 November 2020 c Author(s) 2020. CC BY 4.0 License.

13
Bay of Bengal (Fig. A12). Our simulations show a linear relationship between precipitation over the 335 WEIO and NEA, and precipitation over the Bay of Bengal (Fig. 10c).
We conclude our analysis by investigating the changes in the upper-level (200 hPa) jet (Fig. 11). In the PI experiment, the core of the subtropical jet is located over western Asia and the exit of jet is located over north-eastern China (contour lines in Figure 11). In the MH PMIP simulation, the jet is shifted 340 northwards, with an overall weakening to the south and a strengthening confined to the northward side of the exit of the jet streak (Fig. 11a). The Sahara greening (MH GS ) leads to an accelerated westerly flow at the jet entrance, but an overall slowing down at the jet exit together with a further increase in the northward shift relative to the MH PMIP experiment (Fig. 11b). These changes cause a slight tilt in jet that favours more aloft divergence over northern India and Pakistan as also seen in figure 9b, which in turn 345 favours increased rainfall in the region. The dust reduction (MH RD ) leads to a pattern anomaly very similar to the MH PMIP experiment -albeit weaker (cf. Fig. 11a, c). The effect of the combined forcings (MH GS+RD ) is dominated by the MH GS pattern (Fig. 11d) and in this case the anomalies are even larger than in the MH GS case. This is likely due to the increase in temperature gradient between low and high latitude relative to the MH GS case (not shown).

Model -Proxy intercomparison
To evaluate the model performance when accounting for Sahara greening and reduction in airborne dust concentrations, we compare our simulations to the available marine and terrestrial paleoclimate archives. We focus on the most apparent dissimilarities between the sensitivity experiments and the 355 standard MH simulation (MH PMIP ) where only orbital forcing is considered. While our simulations are centred at 6,000 years BP, they should be seen as indicative of the wet early -middle Holocene rather than a snapshot of exactly 6,000 years BP, which appears to be a period of transition in particular for Indian terrestrial records (e.g., Prasad et al., 1997). 360 https://doi.org/10.5194/cp-2020-142 Preprint. Discussion started: 21 November 2020 c Author(s) 2020. CC BY 4.0 License. 14 Notable differences in summer precipitation between the four simulations occur over western India (Fig. 2), which shows substantially wetter conditions in MH GS and MH GS+RD , compared to the MH PMIP experiment in that region. Nal Sarovar, a brackish lake bordering the Thar desert, appears to have been wetter than today around 6,200 years BP, with a drying tendency towards the end of the MH (Prasad et al., 1997). There is evidence for a substantial pluvial between ~9 and 6 ka farther north in the core of 365 the Thar (Deotare et al., 2004;Gill et al., 2015; and references therein), and a reduced dimension analysis suggests that reconstructed tropical Pacific SSTs alone could have driven a 60% increase in precipitation there during the early Holocene (see figure 5 in (Gill et al., 2017)). However, Gill et al. (2017) inferred winds and the precipitation over India using exclusively a proxy-based reconstruction of the tropical Pacific SSTs, assuming modern teleconnections. The MH PMIP experiment simulates a 370 localized rainfall increase in the region of the Thar desert above 40 -50%, whereas the MH GS+RD suggests a more intense and widespread increase in precipitation (Fig. 2) over the western and northwestern India, even though the monsoonal flow is weaker compared to MH PMIP (Fig. A4). This suggests that the modern teleconnections may not precisely hold in the past and the inferred changes based on only tropical Pacific SST patterns may underestimate the total rainfall changes during the early and 375 middle Holocene over north-western India.
Another region where our simulations show divergent results is south-western coastal India. There, the MH GS+RD experiment shows drier conditions relative to PI, while the MH PMIP shows wetter conditions (Fig. 2). Paleoclimate archives from Nilgiri Hills, in the Western Ghats at the eastern edge of the 380 simulated dry anomaly, suggest that that region was wetter between 12,000 and 10,000 years ago and then during the middle Holocene gradually became drier relative to today (Sukumar et al., 1993).
Hence, accounting for the greening of the Sahara may improve the precipitation anomaly pattern seen during the middle Holocene over India; however, a systematic model validation is not currently possible due to the paucity of available paleoclimate archives and their large uncertainties. 385 With respect to changes in SSTs, the MH GS+RD experiment simulates a warming of the Arabian Sea and Bay of Bengal summer SSTs relative to PI, while little change or a slight cooling is simulated in the https://doi.org/10.5194/cp-2020-142 Preprint. Discussion started: 21 November 2020 c Author(s) 2020. CC BY 4.0 License. 15 MH PMIP experiment (Fig. 4). In the Arabian Sea, proxy evidence for widespread mid-Holocene SST warming is lacking. Dahl and Oppo (2006) showed that the early Holocene (at around 8,000 years BP) 390 was 1.4 ± 1.3 °C cooler than the late Holocene, on the basis of Mg/Ca in the planktic foraminifer Globigerinoides ruber from 12 cores spanning much of the basin. Their only core showing a slight warming at 8,000 years BP (+1.0 °C) was situated off the Horn of Africa, a region with robust warming in all four mid-Holocene simulations. Other G. ruber and Trilobatus sacculifer Mg/Ca records corroborate modest cooling (~0 to -1 °C) during the middle Holocene (around 6,000 years BP) in the 395 eastern Arabian Sea (Anand et al., 2008;Banakar et al., 2010;Govil and Naidu, 2010), with slight warming off the coast of southwest India (Saraswat et al., 2013) and again off the Horn of Africa (Anand et al., 2008). Alkenones document 0 to 1 °C cooling in the northern Arabian Sea during this time period (Böll et al., 2015;Schulte and Müller, 2001), with negligible change off the Arabian Peninsula (Huguet et al., 2006;Rostek et al., 1997) and southwest India (Sonzogni et al., 1998). 400 Regional Mg/Ca and alkenone compilations by Gaye et al. (2018) suggest that no sector of the Arabian Sea was warmer at 6 ka, with the possible exception of south of India, which also warms slightly in all four simulations. Alkenones from the northern Bay of Bengal (Lauterbach et al., 2020) and G. ruber South of the equator west of Sumatra, the MH PMIP simulation produces a stronger cooling (>1 °C) that disappears in the MH GS+RD experiment (Fig. 4). Here alkenones are more consistent with MH PMIP , albeit with a more modest cooling of 0.5 to 1 °C (Li et al., 2016;Lückge et al., 2009). However, G. ruber Mg/Ca indicates negligible change, more consistent with MH GS+RD (Mohtadi et al., 2010). Seasonal 415 differences in proxy carrier production may explain such differences, with Mg/Ca perhaps being more https://doi.org/10.5194/cp-2020-142 Preprint. Discussion started: 21 November 2020 c Author(s) 2020. CC BY 4.0 License. appropriate for comparison to JJAS simulations, as suggested for the equatorial Pacific (Gill et al., 2016;Timmermann et al., 2014). Finally, south of the equator on the western side of the Indian basin off the coast of Tanzania, Mg/Ca reconstruction suggests that SSTs were about 1 to 1.5 °C warmer during the middle Holocene compared to late Holocene (Kuhnert et al., 2014), which is more consistent 420 with the MH GS+RD experiment (Fig. 4d). This record also shows a rapid SST cooling concomitant with an abrupt retreat of the SAM as suggested by a recently published paleoclimate archives from western Yunnan Plateau in southwestern China  and northern Laos (Griffiths et al., 2020).
Such changes are also synchronous with the end of the African Humid Period (e.g., deMenocal et al., 2000), hence our simulations suggest that the changes in vegetation over the Sahara and in airborne dust 425 emissions may have played a key role in shaping the evolution of the SAM.

Discussion and Conclusions
The mid-Holocene was characterised by a strengthening of the northern hemisphere monsoon system (e.g., Sun et al., 2019) due to increased boreal summer insolation. The consequent increase in rainfall led to a greening of several semi-arid and arid regions in Northern Africa and Asia (e.g., Campo et al., 430 1982;Dallmeyer et al., 2013;Fleitmann et al., 2003;Lézine et al., 2011;Tierney et al., 2017), and to a marked reduction in airborne dust emissions (deMenocal et al., 2000;McGee et al., 2013). The largest dust emission decreases are thought to have occurred in Northern Africa, where large tracts of what is today the Sahara Desert were vegetated. Understanding this complex set of interrelated changes can provide insights into the mechanisms of monsoonal variability, and contribute to strengthening our 435 physical understanding of monsoonal changes in climate projections. However, many modelling efforts for the mid-Holocene have focused only on the impact of solar insolation changes as this has been the common protocol for climate simulations of this period (Otto-Bliesner et al., 2016;Taylor et al., 2009Taylor et al., , 2012, neglecting the feedbacks induced by the altered vegetation, soil properties, and associated dust emissions. 440 Indeed, the role of reduced dust emissions during the mid-Holocene on local and global climate has only recently been addressed (Hopcroft and Valdes, 2019;Pausata et al., 2016Pausata et al., , 2017aPausata et al., , 2017bPiao et https://doi.org/10.5194/cp-2020-142 Preprint. Discussion started: 21 November 2020 c Author(s) 2020. CC BY 4.0 License. 17 al., 2020;Sun et al., 2019;Thompson et al., 2019) and it has been shown that airborne dust may play an important role in modulating the intensity and geographical extent of the West African Monsoon 445 (Pausata et al., 2016;Thompson et al., 2019) as well as impacting climate far afield. However, the role of Saharan dust changes in affecting the South Asia Monsoon (SAM) system has not hitherto been investigated. The key goal of the present study is to fill this knowledge gap, by outlining the remote response of the mid-Holocene SAM system to the Sahara greening and associated reduction in airborne dust concentrations. 450 We analyse a set of simulations where the land cover is changed from desert to shrubland over a large part of North Africa and dust concentration over the region is reduced by up to 80% compared to the pre-industrial period (Gaetani et al., 2017;Pausata et al., 2016). We find that a vegetated Sahara -albeit weakening the low-level southwesterly winds -enhances the moisture flux from the Arabian Sea to the 455 northern Indian subcontinent and increases the precipitation in this region compared to a simulation in which only the orbital forcing is considered (Figs. 2, A7 (Figs. 2 and A3). Overall, the SAM rainfall in the MH GS+RD is significantly increased compared to the PI climate as well as to the orbital forcing only simulation (MH PMIP ). The monsoon season is also extended by several months, particularly in the withdrawal phase (Fig. 3). 465 Sun et al. [2019] showed that the greening of the Sahara and a reduction in dust emissions significantly influence the Northern Hemisphere land monsoon precipitation, but the largest impact is on the WAM.
Here, we show that the SAM is significantly affected by both vegetation changes in northern Africa and dust reduction and the remote response is about half of the rainfall change simulated locally over 470 northern Africa (cf. Fig. 2 here with Fig. 2 in Pausata et al. [2016]).

18
A comparison of our simulations with paleoclimate archives points to potential improvements in simulating rainfall over India when including the greening of the Sahara and dust reduction relative to the orbital forcing-only simulation. In particular, our simulations suggest that the vegetation and dust 475 emission changes may have played an important role in affecting the Indian Ocean temperature and shaping the evolution of the SAM during the termination of the African Humid Period. However, no robust conclusions can be drawn in this respect due to the relative paucity of geographically and temporally referenced, quantitative paleo-precipitation data in the region. A similar difficulty is encountered in evaluating the modelled SST changes. Only some paleo-archives point to closer 480 agreement with the MH GS+RD simulation, however, in general the amplitudes of SST changes are small relative to proxy uncertainties, making it difficult to provide a systematic model validation.
Finally, in our experiments we only consider changes in vegetation over northern Africa and its remote impact on SAM. However, proxy archives from the mid-Holocene point to widespread vegetation 485 changes across the globe, with expanded forest cover in Eurasia Tarasov et al., 1998) and greener southern and eastern Asia (Dykoski et al., 2005;Fleitmann et al., 2003;Thompson et al., 1997;Zhang et al., 2014). Swann et al. [2012Swann et al. [ , 2014 show that in their model the remote forcing from expanded forest cover in Eurasia during the mid-Holocene shifts the intertropical convergence zone northward, resulting in an enhancement of precipitation over northern Africa that is greater than 490 that resulting from orbital forcing and local vegetation alone. Using idealized deforestation experiments in the tropics and temperate regions, Devaraju et al. (2015) showed that the monsoonal precipitation changes can be more sensitive to remote than local changes in vegetation. Hence, it is possible that the rainfall changes seen in our study may be further modulated by vegetation changes in Europe and Asia. Therefore, it is critical that the Earth system modelling community conducts a concerted effort to 495 include reconstructed vegetation distributions and dust concentrations when simulating the mid-Holocene climate.

Model Validation
In order to evaluate the performance of the EC-Earth model used here in reproducing the SAM 500 dynamics, we compare our pre-industrial simulation (PI) to temperature and precipitation data from ECMWF's ERA5 reanalysis product, (Hersbach et al., 2020) and gridded observational products. Longterm precipitation rates from ERA5 compare favourably with NASA's TRMM Multi-satellite Precipitation Analysis (Hersbach et al., 2020;Huffman et al., 2010), and over the Indian subcontinent differences between ERA5 and the Global Precipitation Climatology Project (GPCP) gridded 505 observational dataset (Adler et al., 2018) are mostly below 0.5 mm/day (Figs. A1 and A2). Good agreement is also found between ERA5 temperatures and the Climatic Research Unit (CRU) data set (Harris et al., 2020) and ERA5 improves in this respect over previous datasets (Hersbach et al., 2020).
EC-Earth's PI simulation in general underestimates rainfall over the north-eastern Indian subcontinent, and overestimates it over the western side. The model further presents a large cold bias (Fig. A1).

Competing interests
The authors declare that they have no conflict of interest.  Only differences significant at the 95% confidence level using the Student t test are shaded.