Changes in landscape and climate in Mexico and Texas reveal small effects on migratory habitat of monarch butterflies (Danaus plexippus)

The decline of the iconic monarch butterfly (Danaus plexippus) in North America has motivated research on the impacts of land use and land cover (LULC) change and climate variability on monarch habitat and population dynamics. We investigated spring and fall trends in LULC, milkweed and nectar resources over a 20-year period, and ~ 30 years of climate variables in Mexico and Texas, U.S. This region supports spring breeding, and spring and fall migration during the annual life cycle of the monarch. We estimated a − 2.9% decline in milkweed in Texas, but little to no change in Mexico. Fall and spring nectar resources declined < 1% in both study extents. Vegetation greenness increased in the fall and spring in Mexico while the other climate variables did not change in both Mexico and Texas. Monarch habitat in Mexico and Texas appears relatively more intact than in the midwestern, agricultural landscapes of the U.S. Given the relatively modest observed changes in nectar and milkweed, the relatively stable climate conditions, and increased vegetation greenness in Mexico, it seems unlikely that habitat loss (quantity or quality) in Mexico and Texas has caused large declines in population size or survival during migration.

clusters of monarchs 9 ), is driven by linked spatio-temporal generations [10][11][12] .For example, climate conditions at various stages of the life cycle influence seasonal, annual, and possibly longer-term, trends [13][14][15][16] .Weather effects on spring recruitment (primarily in Texas) can influence overall annual population dynamics 13,17 while drought during the fall migration affects monarch lipid levels crucial for overwintering survival 6 , subsequent spring breeding 18,19 , and may influence overwintering population size 20 .
The eastern population declined from the 1990s to 2014, with no statistically detectable trend from 2014 to 2018 21 and recent analyses indicate the population has been relatively stationary, but low, for the last 10 years 22 .Population declines in the mid-1990s and early 2000s were associated with the rapid adoption of herbicidally tolerant GMO corn and soybeans across the midwestern US from ~ 1994 to ~ 2006 and glyphosate application,

Land cover trends
From 2001 to 2020, shrublands dominated landcover in both Texas and Mexico (Table 1, Fig. 2) while crops and developed lands increased through time (Table 1).Notably, developed lands, or most of the developed subclasses, increased >~ 50% (Table 1).
In Mexico, the growth of developed areas and croplands caused a − 2.8% decline in shrubland (Fig. 2a).Forest to shrubland, croplands to shrubland, and barren to shrubland transitions also occurred but at lower rates than LULC away from shrubland.Transition to water as lakes filled and the loss of wetlands to croplands caused a decline in wetland shrub of − 11.9%.
To assure land cover data used in Mexico (Global Land Cover and Discovery (GLAD)) and Texas (National Land Cover Database (NLCD)), measured LULC similarly, we compared statewide trends in gross land cover change across Texas, where both datasets were available from 2000 and 2020 for GLAD and 2001 and 2019 for NLCD.Overall trends were similar (Supplementary results, Fig. S1).County-level differences in levels of gross change could be attributed to differences in the classification systems between GLAD and NLCD.GLAD combines all the NLCD forest classes into a single class, and the GLAD shrubland class included grasslands, pastures, low-intensity and urban open space, barren lands and low-percent-cover forests based on a crosswalk of GLAD and NLCD in Texas (Supplementary Data 1).

Milkweed and nectar resources
Data describing nectar and milkweed resources differed between Mexico and Texas, requiring elicitation of expert opinion for nectar and milkweed resources in Mexico and synthesis of published and unpublished literature in Texas.Across both Mexico and Texas, information on nectar resources was somewhat similar, varying from

Quantitative change in milkweed and nectar resources
Combining LULC with estimates of nectar and milkweed values, resulted in an estimate of a − 2.9% decline in milkweed density in Texas between 2001 and 2019 (starting from ~ 1.04 billion plants in 2001), and no change, or a slight increase, in the milkweed index in Mexico from 2000 to 2020 (Table 3).Nectar resources declined < 0.2% in the fall and spring in both Mexico and Texas.

Geographic patterns of milkweed and nectar change
Most counties (Texas) and municipalities (Mexico) experienced small changes in nectar and milkweed resources (Fig. 3a-c) but there were cases of positive and negative change.Spatial patterns in habitat change indicated more change in Texas than Mexico, and more change (primarily loss) in milkweed than nectar.In Texas, counties along the I-35 highway corridor connecting San Antonio, Austin, and Dallas-Fort Worth, have growing urban areas, and showed declines in spring and fall nectar of 1-5% and milkweed of 4-25%.Counties in central Texas near Abilene exhibited declines in nectar of 1-3% primarily in the spring.One to 5% gains in spring nectar occurred in counties in eastern Texas, while very few counties gained fall nectar.Milkweed density declined 1-5% across most of the Texas study extent, except in the western (mainly southwestern) counties, where relatively large gains of 1-21% occurred, possibly caused by large wildfires converting shrublands to grasslands.Counties showing declines in milkweed density > 15% were along the I-35 highway corridor, the far eastern border near Beaumont, and in a few counties in the north-western parts of the study extent.
In Mexico, eastern municipalities along the US-Mexico border near Reynosa and Matamoros and a small number of additional municipalities demonstrated declines in spring and fall nectar resources.Municipalities www.nature.com/scientificreports/around the overwintering locations had both increases and decreases in spring and fall nectar resources.Milkweed resource value increased 1-5% in a few southern municipalities around the overwintering location, while more northern municipalities, especially those near the US-Mexico border had declines between − 1 and − 7.9%.

Discussion
Changing land cover and climate are fundamental processes driving the environmental conditions affecting species distributions and their population dynamics 54,55 .We measured only small changes in land cover, which resulted in small to no expected changes in nectar and milkweed resources across Mexico and Texas.We also observed few substantive temporal trends in annual climate, aside from an increasing trend in NDVI in Mexico in fall and spring.

Land cover change and monarch habitat
During ~ 2000-2020, the patterns of LULC across the study extent were typical of those observed in other parts of North America where anthropogenic land cover increased, while natural land cover types declined 56 .In Mexico, shrublands declined by ~ 3% while in Texas, none of the natural vegetation classes declined more than 5%, though forests in eastern Texas, where rotational harvest occurs 57 , had the largest declines.The estimated change in milkweed and nectar resources were consistent with the relatively low levels of loss in natural land cover classes and patterns of LULC dynamics affecting vegetation.Our results indicate small, to no, net changes in nectar and milkweed in Texas and Mexico from 1992 to 2021.Milkweed in Texas declined ~ 2.9%, and both milkweed and nectar resources, in fall and spring, showed spatial heterogeneity in trends through time.
In the midwestern US, trends in milkweed were estimated using a combination of field surveys and geographic extrapolations of adoption rates of glyphosate-tolerant corn and soybeans 23,58,59 .Pleasants et al. 58 calculated a loss of 850 million milkweed stems in corn and soybeans, and an additional 11 million lost from changes to land cover in the Midwest, while Flockhart et al. 10 estimated 1.49 billion stems were lost across the entire eastern US monarch breeding range.In non-agricultural lands of the US Midwest, milkweed abundance has likely been stable 59 .
Unlike the US Midwest, we did not detect a change in milkweed in Mexico and estimated a − 2.9% decline in Texas from 2001 to 2019, representing an estimated ~ 29 million plants lost relative to ~ 1.04 billion existing in 2001 across ~ 53 million ha.Our milkweed estimates in Texas are low compared to the 233 million to 1.3 billion plants estimated on just 2.8 million ha of rangelands in Texas 52 .Many of the LC classes in Texas support far fewer milkweeds than rangelands, so we expect our estimates, on a per area basis, to be smaller than Spaeth et al. 52 .Flockhart et al. 10 estimated ~ 1.7 billion stems across the entire southern US, including Texas, indicating our estimate falls between Flockhart et al. 10 and Spaeth et al. 52 .
The impacts of LULC on the eastern migratory population depends on how monarch densities vary across space relative to patterns of increasing or decreasing habitat.Our work indicates there may be geographic concordance between areas in Texas most suitable for monarchs and higher rates of habitat loss.For example, central Texas counties dominated by urbanization and eastern Texas showed larger declines in nectar and milkweed.These same regions may have the highest spring suitability for adults 3 and the urban counties around San Antonio-Austin and Dallas-Fort Worth were most utilized by caterpillars 7 .More formal analysis of these patterns is needed, but this concordance suggests strategies to enhance habitat in urban areas 44,60 and rights of ways in central Texas may help counter the conversion of rangelands and shrublands to urban in these areas.
Mexico showed less overall change in monarch habitat and less spatial heterogeneity in change than Texas, particularly for milkweed.Annual rates of LULC were generally similar across land cover classes for both Mexico and Texas (Table 1).However, in Mexico, experts estimated smaller differences in nectar and milkweed resource values across land cover types than in the US so the largest land cover transitions (shrub to urban, shrub to agriculture) represented relatively smaller changes in monarch resource values in Mexico than in the US.Finally, GLAD in Mexico combines the NLCD grass and shrub classes into a single land cover category and could not estimate changes in habitat caused by transitions between these two cover types; thus, transitions between these cover classes and the difference in value they represent to monarchs could not be captured.

Temporal changes in climate and NDVI
Climate variables showed high levels of annual variability, but the only substantial temporal trend we observed was increased NDVI in both fall and spring in Mexico.The lack of statistically significant trends in Texas runs counter to statewide warming and increased precipitation in eastern Texas observed in longer time series 61 .Our results indicate that during the ~ 30-year time span matching declines in monarch abundance, the climate variables we measured show few to no changes relative to high levels of annual variability.Saunders et al. 20 found a positive relationship between overwintering colony size and fall NDVI in a region covering northeast Texas, Oklahoma, Kansas, Arkansas, Louisiana, and Mississippi.Zylstra et al. 16 , however, focusing on a larger region of the US (90° W to 105° W, 30° N to 40° N) found no relationship between NDVI and overwintering colony size.Both studies used MODIS data, which became available in 2000.We used 30-m LANDSAT derived NDVI from Google Earth Engine (GEE), which is available annually from 1984.Future modelling approaches that include climate variables from Mexico during fall and spring migration would be helpful to determine if they influence population dynamics.
How do our results relate to hypotheses and mechanisms of monarch decline and population dynamics?First, milkweed limitation in Texas does not appear to be a driver of population dynamics, likely because Texas is not heavily dominated by agriculture.Our results, in addition to Spaeth et al. 's 52 estimates, indicate milkweed abundance in Texas has declined, but not to the extent observed in the US Midwest, and that milkweed limitation is unlikely for the first generation of breeding monarchs.
Second, our results indicate loss of migratory habitat in Texas and Mexico may not be a mechanism for potential declines in migration success.We saw little to no signal in the nectar resources that would indicate a decline in fall migration survival and an unexpected increase in NDVI during the fall and spring in Mexico, suggesting the possibility that nectar resources have increased over the last 23 years rather than diminished.
Phenological shifts in milkweed and floral availability may have generated a mismatch between the timing of spring and fall migration 62 and resource availability for monarchs.A separate study is needed to understand if milkweed and nectar phenology has shifted away from the timing of the monarch migration through Mexico and Texas.However, our results indicate that milkweed and floral availability has been predictable for fall and spring migrating monarchs through Mexico and Texas over the last 20-30 years.Migrating monarchs may opt to become residents and take advantage of tropical milkweeds in this region 37 .The monitoring of resident monarch populations in this region will help to test this possibility.www.nature.com/scientificreports/One implication of the migration survival hypothesis is that declines in fall migration survival, large enough to cause a decline in the wintering population, require a commensurate increase in reproductive output in the spring to generate a stable summer breeding population.Furthermore, as the wintering population declined, the difference between the summer breeding population and the winter population increased, so the magnitude of reproductive output necessary to overcome declining fall survival must also have increased year-over-year.
Perhaps the positive trend in spring NDVI in Mexico supports enhanced female spring migration survival from overwintering areas to first-stage breeding grounds and allows greater egg production by females.We are skeptical of this possibility because our data do not show that floral availability has declined in fall and increased in spring, but instead likely increased in both seasonal environments.Ultimately, it is difficult to reconcile how the trends (or the lack thereof) in land and climate change we observed would simultaneously lead to decreased fall migration success yet support or enhance spring migration and reproductive output in equal but increasing measure each year.

Study limitations
The datasets we used can measure shifts in the dominant LULC classes, but the data used cannot adequately map more subtle changes in landscapes like crop margins 48 , nor changes within each land cover class.For example, degraded shrubland sites had different temporal patterns of NDVI than intact sites 63 .Given these issues the LULC, nectar, and milkweed resources analyses were designed to identify large changes in habitat and map them at coarse spatial resolutions 64 .Given the categorical nature of the LULC analyses, these efforts focused on changes in habitat quantity, and we used the climate variables and NDVI, to assess possible changes in habitat quality.
Given the conservation efforts related to butterflies, pollinators, and arthropods in general, developing empirical relationships between host plants, nectar resources, and remotely sensed data would allow broad-scale and historical analyses of these resources and link them to species of interest/management concern.Our Monte Carlo based approach provides a robust framework for integrating LULC data and field collected survey information.Previous studies indicate field campaigns over large regions are feasible 52,59,65 and analyses of the Integrated Monarch Monitoring Program 66 data will hopefully generate milkweed densities for LULC classes over a much broader region than our Texas analysis.
Our results indicate monarch habitat in Mexico and Texas is more intact than the agricultural landscapes of the midwestern US.However, milkweed in Texas has declined, and may continue to do so if existing LULC dynamics persist.As such, efforts to conserve and restore milkweed in Texas, especially in central Texas, where both urbanization and higher densities of monarchs may overlap, may offset some losses.
To the best of our knowledge, this research is the first monarch study to analyze trends in climate and NDVI across migratory habitat in Mexico and the observed increasing trends in NDVI were unexpected.Incorporating this information into models of eastern monarch population dynamics would help us better understand if and how spring and fall conditions in Mexico affect population dynamics.A next step is to understand how much habitat is currently protected from future land cover changes, resilient to climate change, and if these protected areas would be sufficient to support monarch migration and reproduction and generate similar benefits to other species.We anticipate continuing changes to climate, so understanding climate-related changes in monarch, floral, and milkweed distributions will inform decisions regarding how and where to best conserve existing migratory habitat.

Materials and methods
Our methodological choices were based, in part, on data availability.For example, milkweed and nectar densities are not well-sampled across our study extent, and land cover data vary between the US and Mexico 64 .We used land cover datasets from Mexico and Texas to estimate land cover change, then used Monte-Carlo methods to link changes in land cover to changes in monarch habitat (spring and fall nectar and milkweed resources) using empirical data and surveys of experts.In addition, we assessed temporal trends in climate variables and NDVI across the study extent in Mexico and Texas to determine if climate conditions that affect plant growth and nectar production, as well as an index of vegetation condition, changed through time.NDVI measured vegetation greenness and recent studies suggest it, or the related Enhanced Vegetation Index (EVI) correlates with nectar availability 65 .Winter and spring NDVI values have also been shown to influence butterfly population dynamics, such as painted lady butterfly (Vanessa cardui) migrations in Africa 67,68 .

Study time frame and extent
We analyzed LULC data from ~ 2001-2022 as well as climate and NDVI data from 1992-2021 given differences in when these data became available.We provide specific time periods and data used for Mexico and Texas below.We demarcated our study extent in Mexico by intersecting the extent from 4 with the municipalities in Mexico (Fig. 1).For Texas, we used counties in the US Fish and Wildlife Service (USFWS) Monarch Butterfly Conservation Units 69 but added coastal counties in Texas (see Ref. 70 and locations of monarchs observed in the Journey North data (https:// journ eynor th.org/)).This extent generally agrees with the study extent in Texas used by Ref. 7 and their modelled kernel densities of monarch larvae.We combined the municipalities in Mexico and counties in Texas for a final study extent (Fig. 1).

Land cover data
We previously reviewed land cover datasets that could be used in comparisons across Mexico and Texas 64 .For this study, we wanted land cover data with high levels of spatial resolution to track change in habitat that also covered the time series of overwintering abundance in Mexico (1994-present).This requirement necessitated using datasets unique to Mexico and Texas.For Mexico, we used GLAD data for 2000 and 2020 71 to track LULC patterns in 8 land cover classes developed from the GLAD land cover strata.GLAD data are currently available for these two years at 30-m resolution.Prior to selecting GLAD, we performed LULC analyses using the European Space Agency Land Cover (ESA-LC) product (1992-2018).This 22-class raster dataset has a 300-m resolution and is produced annually, which seemed ideal for our purposes.However, our analyses indicated the 300-m resolution data did not detect land cover change at acceptable levels of precision (Supplementary results, Fig. S4).We also considered the 30-m resolution North American Land Change Monitoring System (NALCMS) data but at the time we initiated this study, they were only available in 2010 and 2015, though data for 2020 recently became available.Relative to GLAD, both NALCMS and ESA-LC overestimated agriculture in many parts of our Mexico study extent and GLAD maps best-matched aerial imagery (Supplementary results, Fig. S4).
Finally, we compared GLAD and NLCD in Texas, where both datasets overlapped, to check how similarly they detected land cover change.For each dataset, we developed maps of percent total gross change in land cover for the bookend years across counties (GLAD, 2000-2020, NLCD, 2001-2019).

Estimating land cover change
For both Mexico and Texas, we calculated all transitions from one land cover class to another, including no change, by comparing raster data layers in 2001 and 2019 in Texas and 2000 and 2020 in Mexico.For each time interval, we calculated gains, losses, net change, gross change, and annual change for each LULC class.Net change is gain minus loss whereas gross change is gain plus loss.Annual is annualized net change.We summarized land cover change for each municipality (Mexico) or county (Texas) and mapped this change.In addition, we generated chord diagrams of land cover transitions, to visualize the dominant 'many to one' changes in land cover classes from the first to last year of the time series.

Linking land cover change to monarch habitat
Our Monte-Carlo simulations required statistical distributions of milkweed and nectar resources for each land cover class.In Texas, we used the national spring and fall floral resource values from Ref. 72 to translate land cover to floral resources.Koh et al. 72 surveyed experts to estimate floral resource values for each season, NLCD class, and crop type in the Cropland Data Layer (CDL) (See SI. pnas.1517685113.sd01.xlsx).NLCD only has a single crop type, so we calculated a weighted average of Koh's floral values for crops, where weights were the average percent for each CDL crop class making up the NLCD crop class across each of the main NLCD years (2008,  2011, 2013, 2016, Supplementary Data 2).We excluded the grass category due to issues with parsing out sod/seed/ grass seed from NLCD grass/hay categories.This arrangement allowed us to use the single cropland class in the NLCD while assigning a nectar value that represented the dominant crops grown across the Texas study extent.
In Mexico, we lacked empirical data of milkweed and nectar abundance among land cover classes.As an alternative, we duplicated the survey developed by Koh et al. 72 to elicit expert opinion regarding fall milkweed and spring and fall nectar resources in Mexico.The survey (Supplementary results, Figs.S5-S7) was designed and administered by coauthor FB via SurveyMonkey to 30 monarch butterfly experts who participated in the Red Nacional de Monitoreo de la Mariposa Monarca en México National Commission of Protected Area (Working Group for the Conservation and Monitoring of the Monarch Butterfly Flyway).This monitoring program has generated ~ 45,000 monarch sightings, and records monarch behavior on flowering plants, the overnight resting sites on the migratory route, the number of butterflies, and other variables 8 .These biologists were very familiar with INEGI, Mexico's primary land use and vegetation data 73 , https:// www.inegi.org.mx/ temas/ usosu elo/).The INEGI data contain ~ 135 vegetation and land cover classes, which was too complex for the survey.Instead, we used the land cover classes from a 14-class generalization of the more detailed INEGI vegetation classification system 74 .
The survey asked respondents to provide information about the quality of nectar resources in the spring and the fall, as well as the density of milkweed plants in the fall, for each of the 14 generalized INEGI land cover classes.For nectar, we followed the wording of questions used by Koh et al. 72 , and respondents scored each land cover class with 5 levels of relative resource availability labelled as "0.10", "0.25", "0.5", "0.75", "0.95".Even though the 5 levels have quantitative labels, the question used a Likert, or rating scale and the 5 categories represent experts' subjective opinion about the quality of the nectar resource available for monarchs.The question does not represent a density estimate such as nectar production per unit area.For milkweed, we modified the nectar question slightly, to better address density, because we had density information in Texas.We initially considered asking experts to provide a specific density estimate (plants /ha) for each land cover class to match the units of Texas data, but some experts did not feel comfortable providing this level of detailed information.Instead, we asked respondents to categorize the density of milkweed, as percent cover, across the 5 categories, where each category label represented a density class (0.10 = 10%, Fig. S8, Supplementary Results).Given this limitation, the milkweed density data in Mexico provide relative differences across land cover classes, but do not allow us to measure plants per ha, or convert density to an estimate of total plants, as we did in Texas.
The data were normally distributed so means and standard deviations were calculated across respondents for each land cover class.Anonymized survey results are available in Supplementary Data 3.
Different land cover datasets are rarely fully in agreement.To incorporate the uncertainty between the 14-class INEGI used in the survey and the 8-class GLAD used for land cover change, we used agreement matrices between INEGI and GLAD to translate the survey-based scores for milkweed and floral resources, based on INEGI classes, to the GLAD classes (for an example with the GLAD forest class, see Supplementary results, Table S1).
Vol:.( 1234567890 For spring milkweed in Texas, we analyzed milkweed density (plants/ha) data (n = 613) from a variety of sources compiled by the Xerces Society for U.S. Geological Survey (these data are in Supplementary Data 5), 33 of the 613 samples were proprietary IMMP data collected in Texas (provided by the Monarch Joint Venture and USFWS).Xerces contacted 140 state biologists and university researchers to generate the data.These data represent different field efforts designed to sample milkweed density.They were collected in 2004-7, 2009, 2011,  2015, and 2016.Data were filtered to include samples that were in Texas, had geographic coordinates and were density estimates (plants/ha) not presence/absence data.In addition, we excluded data with comments that called their quality into question.For example, some excluded records noted "student maxed out at 15" indicating a plot had more than 15 plants, but the field technician stopped counting at 15.The remaining 613 samples were overlaid with NLCD, then manually checked with aerial imagery.The NLCD landcover classes were adjusted when aerial imagery indicated a misclassification of the land cover.This adjustment happened mainly along roads, where natural lands were misclassified as developed.
Across the entire state, and within each land cover class, the milkweed density data were highly skewed, dominated by zeroes (60-91% across land cover classes).We selected the best fitting distribution for each land cover class, by comparing zero-inflated binomial, exponential, negative binomial, and zero-inflated Poisson distributions (Supplementary Data 7).Despite > 600 samples, six land cover classes had low sample sizes for estimating milkweed densities.Five of the classes, Open Water, Barren Land, Cultivated Crops, Woody Wetlands and Herbaceous Wetlands, had no milkweeds present in the samples.Lacking better data, we set these values to zero.The final land cover class with sparse data, Developed Medium Intensity, was estimated as the midpoint between the Low-and High-Density Developed classes.
Note that the survey-based estimate of floral resources (spring and fall in Texas and Mexico) and milkweed (spring, Mexico) are unitless, representing expert opinion of suitability on a scale from 0 to 1. Adding these across pixels creates total values that are difficult to interpret, so we report these data as percent change, and they represent relative differences through time.In Texas, the spring milkweed data represent actual plant density (plants/ha), and we report changes in this value, as total plants, as well as percent change.

Monarch habitat change
To estimate extent-wide changes in monarch habitat in Mexico from 2000 to 2020 and Texas from 2001 to 2019, we combined the area estimates of change in each land cover class with the distributions of nectar value or milkweed abundance for each land cover class using Monte-Carlo simulation.For milkweed density in Texas, we used the best fitting distribution among those we evaluated.For Texas nectar, Koh et al. 72 reported parameter values from continuous beta probability distributions, which we used.In Mexico, as noted above, survey data were normally distributed and we used means and standard deviations for each LC class to generate distributions for spring milkweed and spring/fall nectar.
For 2000 and 2020 in Mexico, and 2001 and 2019 in Texas, we first multiplied the areal amount of change by a random draw from the distribution of floral resource or milkweed value for the respective season and land cover class.We then summed the nectar or milkweed value across all land cover classes to estimate an extent-wide value for each study year and then subtracted the first year from the last year to estimate change from 2000 to 2020 in Mexico and from 2001 to 2019 in Texas.This process was repeated 5,000 times and we report the mean and standard error of the change and the percent change from the 5,000 draws from the start and end years in Mexico and Texas.
To visualize spatial patterns of change across the entire study extent, and identify areas of increasing or decreasing habitat, we mapped the amount of change in nectar and milkweed resource value across municipalities (Mexico) and counties (Texas)

Changes in climate and NDVI
To estimate climate trends in the spring and fall, we gathered data from Google Earth Engine (GEE), summarized the data into seasonal median values for each year, and used linear models to estimate trends through time.We developed Python 75 scripts to extract TerraClimate data and Landsat derived NDVI from 1992 to 2021, with the GEE API and geemap Python library (https:// geemap.org).All data were sourced from the GEE database (https:// devel opers.google.com/ earth-engine/ datas ets).Data from TerraClimate (https:// www.clima tolog ylab.org/ terra clima te.html) included precipitation (precip, mm), maximum temperature (tmax, °C + 10), minimum temperature (tmin, °C + 10), and Palmer Drought Severity Index (pdsi, unitless), at a resolution of 4638.3 m.NDVI data came from (https:// www.usgs.gov/ cente rs/ eros/ scien ce/ usgs-eros-archi ve-veget ation-monit oringnoaa-cdr-ndvi), and provide a measure of vegetation greenness, an important indicator of vegetation health.NDVI data have a 30 m resolution but were summarized at a resolution of 300 m.
We split the climate and NDVI data into fall (September, October, and November) and spring (March, April, and May) and analyzed these separately.For each year, we calculated a median across all available estimates www.nature.com/scientificreports/within a season using the median reducer from GEE.For example, satellites collect NDVI every 8 days, so for a 3-month season (~ 90 days), we had ~ 11, 8-day values from which we calculated a median.This was done for each pixel, for each year, to create a median composite NDVI image across the entire study extent.We repeated this process, for each variable, for every year from 1992 to 2021.
For each year, we summarized the map across the study extents for Mexico and Texas by calculating a median across all pixels.The median value was then assembled into a time series, forming a set of climate and NDVI data that covered ~ 30-year period in Texas and Mexico.We used simple linear regressions to check for temporal trends where the climate variable or NDVI was the dependent variable and year the explanatory variable.
In addition to the supplemental data included in the paper, data and code associated with the analyses are available at Ref. 76 .

Figure 1 .
Figure 1.Study extent in Mexico (with municipalities) and Texas (with counties) created using ArcGIS Pro 3.2.1,Copyright © 1995-2023 Esri (https:// www.esri.com/ en-us/ arcgis/ produ cts/ arcgis-pro/ overv iew).MBBR is the outline of the Monarch Butterfly Biosphere Reserve and represents the general region where monarchs overwinter.

Table 3 .
Net change in milkweed and nectar resources from 2000-2020 in Mexico and 2001-2019 in Texas.Data are presented as means (SE in parenthesis) of both raw values and total percent change relative to the first year of the time interval.Note that total percent change is calculated from the start and end years of the interval for each country/state and is not an annual rate of change.For Texas, milkweed is reported as total plants across the study extent.All other values were derived from expert opinion and represent resource value change per km 2 .

Figure 3 .
Figure 3. Patterns of percent change across the study extent for spring nectar (A), fall nectar (B), and Spring milkweed (C).Resource change was summed at the county level in Texas and the municipality level in Mexico created using ArcGIS Pro 3.2.1,Copyright © 1995-2023 Esri (https:// www.esri.com/ en-us/ arcgis/ produ cts/ arcgis-pro/ overv iew). https://doi.org/10.1038/s41598-024-56693-z https://doi.org/10.1038/s41598-024-56693-zwww.nature.com/scientificreports/ Because the milkweed resource values data were on different scales (field-based measurements in Texas vs expert opinion in Mexico), we fit the estimated mean change of milkweed density in Mexico to the range established in the Texas field data.To do so, we calculated a z-score of the Mexico milkweed density change values, their means, and their standard deviations, and then subsequently de-standardized those values by using the mean and standard deviation of the Texas milkweed density change values.The result is a range of milkweed change values across Texas and Mexico all within the same range of values. https://doi.org/10.1038/s41598-024-56693-z

Table 1 .
2and cover change in Mexico and Texas (km2).Net change is Gain minus Loss whereas Gross is Gain plus Loss.Annual is annualized net change.Land cover class labels follow the datasets used for Mexico and Texas.Herbaceous is shortened to herb; Dev.developed.comparativelylowvalue(~0.2) to moderately valuable (~ 0.5) among land cover classes (Table2).Differences between fall and spring in the value of individual land covers to monarchs were greater in the more temperate US compared to Mexico.Mxico experts considered agriculture, scrublands, and urban as supporting the highest spring nectar resources, with forests and shrubs having the highest values in the fall.For Texas, shrublands, grasslands, and herbaceous wetlands had the highest spring nectar values, while deciduous forests and shrublands were scored highest in the fall.Experts in Mexico gave similar milkweed values (0.23-0.29) across all land cover types except the lower-valued wetland herbaceous land cover.Across Texas, field-collected milkweed density varied an order of magnitude (~ 6 to 60) across land cover classes.Grasslands, deciduous forests, and developed open space had the highest milkweed densities, while evergreen forests, and dense urban development had the lowest (Table2).
Figure 2. Chord diagram showing land cover changes across the study area in Mexico (A) from 2000 to 2020 based on Global Land Cover and Discovery (GLAD) data and Texas (B) from 2001 to 2019 based on National Land Cover Database (NLCD).For Texas, individual Forest, Developed, and Wetland classes were combined into single categories to enhance legibility.The same color leaving a cover class to another represents a transition from that class to the other.GLAD classes were labeled akin to NLCD for ease of interpretation across figures.A crosswalk between GLAD and NLCD is available in Supplementary Data 1.Vol.:(0123456789) Scientific Reports | (2024) 14:6703 | https://doi.org/10.1038/s41598-024-56693-zwww.nature.com/scientificreports/

Table 2 .
Nectar and milkweed resource values by land cover classes for Mexico and Texas, including means (and standard deviations) from the modelled distributions.Values for Mexico were estimated from surveys.Water was zero for both countries.Values for Fall and Spring Nectar in Texas were from Koh et al. 2016.Milkweed values in Texas were estimated from samples as plants/ha.Land Cover Class labels follow the datasets used for Mexico and Texas.

cover class Fall nectar Spring nectar Milkweed
Vol:.(1234567890) Scientific Reports | (2024) 14:6703 | https://doi.org/10.1038/s41598-024-56693-z )Each GLAD class comprised different proportions of the INEGI land cover classes.We used these proportions to calculate a weighted average of both the mean and standard deviation of the survey-based scores for each GLAD class.We averaged the agreement/disagreement proportions between GLAD and INEGI versions (series 3, 4, 5 and 6) and each year of GLAD, and used these values when calculating the weighed means and standard deviations of the survey-based nectar and milkweed scores (Supplementary Data 4).