Impacts of climate change on winter wheat and summer maize dual-cropping system in the North China Plain

Climate change exacerbates the frequency of agricultural-relevant climate extremes, which could threaten crop growth and food production. The North China Plain (NCP), as one of the most important food production bases in China, is dominated by winter wheat and summer maize. The response of this dual-cropping system to climate change has not been thoroughly evaluated. In this study, the remote sensing normalized difference vegetation index (NDVI) was used to extract the dynamic phenology of winter wheat and summer maize and characterize crop growth status. The meteorological data from 1982 to 2015 were used to examine the mean climate factors and four typical climate extremes (including heat stress, spring frost, meteorological drought, and extreme wet events) associated with phenology shifts. Then, the effects of climate changes on winter wheat and summer maize growth were explored by a panel regression model. In the NCP during 1982–2015, the winter wheat growth exhibited no significant sensitivity to the four climate extremes, and only extreme wet event exerted a significant impact on summer maize growth. The insensitivity of crop growth to climate extremes may benefit from widespread irrigation, improved cultivars and agricultural management (e.g., topdressing and insect pests control). In addition, over the last 34 years, mean climate conditions, especially average temperature, solar radiation and vapor pressure deficit, generally made more contributions to the variations of wheat and maize growth than climate extremes, indicating that mean climate conditions dominated crop growth changes in the NCP. Our findings highlight the possible effects of climate change on crop growth of regional dual-cropping system and provide a critical foundation for future effective measures to ensure regional food production.

pollen (Lobell et al 2013), and changing grain number (Liu et al 2014, Barlow et al 2015, Akter and Islam 2017. For example, Lobell et al (2011a) concluded that each degree day above 30°C during the growing season reduced the final maize yield by 1%-1.7% in Africa from 1999 to 2007. Using the statistical model of rice yield responding to temperatures at stations, Shi et al (2015) found that rice yield in South China decreased by 1.5%-9.7% due to the post-heading heat degree days above 35°C from 1981 to 2010. On the other hand, cold stress can affect crops both in vegetation and reproduction periods (Thakur et al 2010). Wheat is one of the most sensitive crops to cold stress, as its cold tolerance drops rapidly after seeds breaking through the soil surface and greening up in spring. Li et al (2015) reported that a five-day spring frost at the jointing stage could lead to 5%-14% yield loss of winter wheat based on a two-year controlled field experiment. The temperature control field experiment on three winter wheat cultivars also pointed out that the more serious damage to kernel weight and yield per plant occurs when the subfreezing temperature drops from −1°C to −9°C during the jointing stage (Wu et al 2014).
Drought is one of the severest extreme weather disasters that pose water threats, damage crop growth processes and constraint crop production (Piao et al 2010, Lesk et al 2016, Wang et al 2018. The influence of droughts on crops includes affecting stomatal conductance, increasing floral abortion, shortening grain filling duration (Prasad et al 2008), and advancing crop harvest date (Ren et al 2019). For example, Lu et al (2020) found that drought stress qualified by the standardized precipitation evapotranspiration index (SPEI) occurring during the reproductive stage had more serious effects on non-irrigated crops, and one-unit decrease in SPEI could cause about 20% more damage to non-irrigated soybeans than irrigated soybeans in the United States from 1950 to 2016. In addition, excessive precipitation can affect crops by delaying planting, damaging young crops, causing anoxia, root diseases and soil erosion (Urban et al 2015, IPCC 2019. According to the historical maize yield and meteorological data in the United States from 1981 to 2016, a comparison of yield in extreme precipitation years with long-term yield trends demonstrated that excessive rainfall could reduce maize yield by up to −34%, comparable to the impact of drought (Li et al 2019b).
Besides climate extremes, research on mean climate factors (e.g., temperature, precipitation, solar radiation and vapor pressure deficit) during crop growing season is indispensable (Schlenker and Roberts 2009, Tao et al 2013, Lobell et al 2014, Tao et al 2016, Tao et al 2017, Li et al 2019a. A residual trend analysis of climate factors and crop vegetation index in West Africa demonstrated that the combined effect of precipitation, solar radiation and land surface temperature could explain about 40% of the crop productivity variation (Mechiche-Alami and Abdi 2020). Under global climate change, rising temperatures may be detrimental to crop yield by accelerating respiration rates, shortening growing period and reducing organic matter accumulation (Lobell and Field 2007, Tao et al 2012, Xiao and Tao 2014, Zhao et al 2017, Chen et al 2018, but they may also bring benefits if temperature is below the threshold temperature level (Tao et al 2013, Zhao et al 2016, Li et al 2019a. For example, the analysis of eddy covariance flux and meteorological observations in a maize field in Northeast China indicated that the increasing temperature shortened the carbon uptake period but improved the maximum net ecosystem productivity and the complementary effects determine cumulated maize productivity (Zhou et al 2021). As a critical water supply, total precipitation during growing season is positively associated with crop production (Lobell and Field 2007, Lobell et al 2011b, Zhang and Yao 2013, while the impact could be weakened in irrigated cropland (Xiao and Tao 2014). Vapor pressure deficit, measuring atmospheric dryness, is another key element influencing crop growth and development. Increased vapor pressure can trigger stomatal closure and limit water use, leading to yield loss Troy 2018, Yuan et al 2019). From 1995 to 2012 in the America Midwest, the regression model between maize yield and vapor pressure deficit demonstrated that the vapor pressure deficit reduced maize yield by −15% to −28% (Lobell et al 2014). Hsiao et al (2019) also showed that the increased vapor pressure deficit associated with 2°C warming could lower maize yield by 12.9%. Solar radiation is another necessary factor affecting plant photosynthesis and biomass accumulation , and its diverse impacts on crop yield were uncovered in previous studies (Chen et al 2013, Tao et al 2015, Tao et al 2016. For example, wheat and maize yields were likely to increase in the southeast but decrease in other regions of China under increasing solar radiation, as high solar radiation can increase evapotranspiration and reduce rain-fed crop yields (Tao et al 2012). Decreased solar radiation by 10% resulted in more than 20% decrease in wheat yield at some stations in North China during 1981(Tao et al 2017.
The North China Plain (NCP) is the primary grain production area of China. Winter wheat and summer maize are the most important crops on the NCP, accounting for 75% and 35% of the national wheat and maize production, respectively (Mo et al 2009. Various measures, such as irrigation facilities, fertilization and crop variety improvement have been implemented to maintain regional wheat and maize production, making it a representative winter wheat and summer maize dual-cropping area (Shen et al 2013). A series of studies have investigated the relationship between climate extremes and crop production in the NCP. Shao et al (2021) found that temperatures exceeding 35°C decreased the grain number at the ear tip and reduced the final maize yield through field experiments on summer maize. Based on the observed data of agricultural meteorological stations from 1981 to 2009, Tao et al (2015) investigated that each unit increase in high temperature (>34°C) degree days (HDD) in the reproductive period reduced wheat yield by −0.17% to −0.34% in most of the NCP. Liu et al (2014) used a similar method to demonstrate that post-heading heat stress (heat degree days above 30°C) exhibited a significant negative effect on winter wheat yield with the sensitivity of yield to heat stress of −0.01. Xiao et al (2018) showed that wheat yield in the NCP declined by about 1.5%∼2.1% for each unit increase in the spring frost index. Using Pearson correlation analysis,  found that the droughts quantified by the 3-month SPEI in April severely affected the simulated winter wheat yield in the central NCP. Li et al (2019a) indicated that drought quantified by the Palmer drought severity index negatively affected summer crops but positively affected spring crops. These previous studies mainly analyzed a single climate extreme (Lobell et al 2012, Lobell et al 2013, Liu et al 2014, Chen et al 2018, Wang et al 2018, Xiao et al 2018. The relative importance of different climate extremes on crop growth is still unclear. Except for climate extremes, average temperature was suggested to have a negative relationship with spring crops (e.g., winter wheat), while a positive relationship with summer crops (e.g., summer maize) (Li et al 2019a). Chen et al (2017) also pointed out that one unit increase in average temperature enhanced winter wheat yield by 2.1% in the northern NCP and decreased yield by 4.0% in the southern NCP. In addition, solar radiation was positively correlated with winter wheat and summer maize yield at most of the agro-meteorological stations in the NCP (Tao et al 2016(Tao et al , 2017. Maize yield was also significantly positively correlated with precipitation, and its sensitivity to precipitation was much lower than that to temperature and solar radiation within our study region according to field observations from stations (Tao et al 2016) Although previous studies have explored the impacts of climate change on crop growth, the relative contribution of climate extremes compared with mean climate factors to the dual cropping system in the NCP remains unknown.
As the main distribution area of the winter wheat and summer maize rotation system in China, the NCP was selected as the study area. Combing the long-term remote sensing and meteorological data, the dynamic phenology of winter wheat and summer maize were identified, and the agricultural-relevant mean and extreme climate conditions during the crop growing season were obtained. The specific objectives of the study were to (1) examine the temporal and spatial changes in the mean and extreme climate factors during 1982-2015, and (2) quantify the impact and contribution of climate factors on the growth of winter wheat and summer maize rotation system in the NCP by using a panel regression model.

Study area
The NCP is located in the east of China (32°N-40°N and 114°E-121°E). With a typical temperate monsoonal climate, the annual mean temperature is 14°C-15°C and precipitation ranges from 500 mm to 1000 mm, mainly occurring from June to September (Mo et al 2009, Wang et al 2018. The rotation of winter wheat and summer maize is the dominant cropping system in the region, which is mainly distributed in the piedmont plain of the Taihang Mountains, the irrigation districts along the Yellow River and the southwest part of the NCP (figure 1). Winter wheat is cultivated between October and June in the next year and the following summer maize grows from July to September. Due to highly intensive towns and built-up areas, the cropland in the NCP is relatively fragmented (Li and Lei 2021). Winter wheat and summer maize are sown and harvested by local farmers using tractors or combine harvesters. New cultivars of wheat and maize were frequently released by local crop seed stations and authorities (Li et al 2019a). Fertilizers and pesticides have also been applied during the growing season to prevent pests and diseases and ensure crop production (Mo et al 2009, Xiao andTao 2014).
In the NCP, irrigation is necessary for grain production (Shen et al 2013. Groundwater and river withdrawal are the main sources of irrigation in the north and the area along the Yellow River (Mo et al 2009, Shen et al 2013, while crops are less irrigated in the southern NCP owing to the relatively sufficient precipitation. In the rotation system, winter wheat consumes over 70% of the irrigation water  and one to three irrigations are required. Since summer maize grows in the rainy season in a year, its water consumption can basically be met by precipitation, needing less than two supplemental irrigations. Multiple irrigations are implemented, the irrigation amount is still insufficient to fully meet the water requirements for crop growth at the entire NCP level (Lei et al 2015).

Data
The data used in this study contained remote sensing images, statistical yield data, observed weather and phenology data, and the winter wheat-summer maize distribution map.
The third-generation Global Inventory Modeling and Mapping Studies (GIMMS3g) normalized difference vegetation index (NDVI) dataset with a spatial resolution of 0.08°and a temporal resolution of 15 day from 1982 to 2015 was obtained from the NASA Ames Ecological Forecasting Lab (https://ecocast.arc.nasa.gov/data/ pub/gimms/). This dataset has been widely used to detect phenological stages, monitor vegetation status and forecast crop yield (Wang et al 2017, Li et al 2019a). In this study, GIMMS3g NDVI was applied to extract critical phenological dates and then to calculate the seasonal average NDVI value for winter wheat and summer maize. To assess the correlation of NDVI values and crop yield, we also collected a series of city-level statistical annual yield data of winter wheat and summer maize in the study area. The cities with long-term statistical yields (no less than 28 years) were selected and were mainly located in the central and southwest of the NCP (figure 1) and the statistical data are available in the provincial statistical yearbook (http://data.cnki.net/).
The daily weather data included maximum, minimum, and mean temperature (°C), precipitation (mm), sunshine hours (h), minimum relative humidity (%) and minimum ground surface temperature (°C) during 1982-2015 and were collected from 363 meteorological stations in and around the NCP provided by the China Meteorological Data Service Center (CMDC) (http://data.cma.cn/en) (figure 1). The phenology records of winter wheat and summer maize, such as sowing, green-up, heading and maturity spanning from 1991 to 2013 were from 118 agro-meteorological stations within the NCP (figure 1). The data were collected from the monitoring field at the station, and the agricultural management was consistent with the local traditional practices to ensure its representativeness. The data could also be found in CMDC.
The distribution map of the main cultivation region of winter wheat-summer maize in the NCP was created using the supervised classification results of the MODIS NDVI images at a 0.0025°resolution during 2001-2018 (Li and Lei 2021). To unify the spatial resolution, the planting area map was resampled to 0.08 degree.

Data pre-processing
To remove the noise caused by cloud contamination or poor atmospheric conditions, the GIMMS3g NDVI time series were reconstructed by the Savitzky-Golay filter of the TIMESAT software (Eklundh and Jönsson 2009). The parameters in the Savitzky-Golay filter including the smoothing window size, number of envelope iterations and spike method were set to two, one and two (seasonal-trend decomposition using loess replacement), respectively. All daily meteorological data of stations were interpolated into 0.08°grid data using the distance and direction weighted mean method (details can be found in Yang et al (2004)). In addition, we calculated the average values of recorded green-up, maturity date of winter wheat, and maturity date of summer maize in 1991-2003 at each agro-meteorological station. Then, the station-based results were transformed to grid values with the same interpolation method as the weather data. The three multi-year average phenology maps provided a reference for extracting key phenological points of winter wheat and summer maize from the remotely sensed NDVI curve. To eliminate the influence of other cropping systems, we applied the resampled winter wheat and summer maize distribution map to select winter wheat and summer maize grids for subsequent analysis.

Identifying winter wheat and summer maize phenology
We discriminated key phenological dates in each planting area grid from 1982 to 2015 based on the NDVI curve, including the start of season (SOS), peak of season (PEAK) and end of season (EOS) of winter wheat and summer maize. For winter wheat, we did not identify the actual sowing date in early autumn but the SOS of the NDVI profile in the spring season, seeing that the rapid growth after this date dominates photosynthetic activity and crop development ( show that winter wheat reaches the heading stage around April to May, its PEAK was identified as NDVI reached the maximum value (NDVI max ) during the 90th day of the year (DOY90) to DOY140 (roughly corresponding to March 31st -May 20th). Similarly, summer maize peaked from DOY191 to DOY 263 (roughly corresponding to July 10th -September 20th) (Xin et al 2003, Lu et al 2014, Li et al 2019a. We also extracted the minimum NDVI point in a year from PEAK of winter wheat to PEAK of summer maize as the distinct valley point between two growing seasons. Third, we successively retrieved the SOS of winter wheat, EOS of maize, and EOS of winter wheat referring to an effective extraction scheme combining the vegetation index curve slope, thresholds, and multi-year average phenology records in Chen et al (2018). The flowchart for identifying the SOS of winter wheat is shown in figure  S1 (available online at stacks.iop.org/ERC/4/075014/mmedia) and details of the extraction scheme are in Text S1. Finally, from the minimum NDVI point (i.e., the valley point) to PEAK of summer maize, the SOS of summer maize was measured as the start day of a continuous positive slope above two thresholds (Araya et al 2018). One threshold equaled NDVI min (NDVI value at the valley point) plus 20% of the difference between NDVI max (NDVI value at PEAK of maize) and NDVI min , the other was set 0.2 which represents bare soil (Araya et al 2018). The growing season length (GSL) of wheat and maize was defined as the duration from SOS to EOS. To evaluate the accuracy of the extracted phenology, we selected the nearest four winter wheat-summer maize grids (within 20 km) around each agro-meteorological station for each year to calculate the average phenological dates. Then, the averages (i.e., SOS, PEAK and EOS for wheat and maize) were compared with the individual observed phenology (e.g., green-up, jointing, heading date of winter wheat and the seventh true leaf, jointing dates of summer maize) at the corresponding station.

Calculation of climate and NDVI variables
In the study, we used the extreme degree days (EDD,°C ·h) and spring frost temperature-drop-rate (TDR,°C ·d −1 ) to quantify extreme heat stress and cold stress, respectively. The standardized precipitation index (SPI) was calculated to denote the meteorological droughts and extreme wet events.
The extreme degree days (EDD) was the accumulation of degree days when the temperature exceeded a threshold (Lobell et al 2012, Lobell et al 2013, Shi et al 2015. The index characterizes both the intensity and duration of heat stress and is effective to evaluate the effects of heat stress on wheat and maize (Lobell et al 2012, Lobell et al 2013, Liu et al 2014, Zaveri and Lobell 2019. The EDD in each winter wheat and summer maize grid during the growth season (i.e., from SOS to EOS) was obtained as follows: where DD t is the degree day for hour t, t is the hourly time step, N is the total number of hours in the growth season. T h represents a high temperature threshold and is set to 25°C for winter wheat and 32°C for summer maize (Deryng et al 2014). T t is the hourly temperature and is calculated by fitting a cosine function to daily maximum and minimum temperatures .
Because we focused on the growing season of winter wheat after green-up, the extreme cold event in winter was not considered. Considering that the temperature-drop-rate can well depict the severity of cold stress, the spring frost temperature-drop-rate (TDR) representing the average drop rate of all the consecutive cooling periods when the minimum ground surface temperature falls below zero during the vegetation phase (from SOS to PEAK) of winter wheat was obtained by equation (3) where GST min is the minimum ground surface temperature and GST min,di is the minimum ground surface temperature on day di. When GST min on day d1 is below 0°C (China Meteorological Administration 2008, Bao et al 2012, Xiao et al 2018) and GST min continues to decrease in the following days till that GST min reaches the lowest value on dk, frost occurs and these days are regarded as one cooling period. The beginning and ending days of this cooling period are d1 and dk, respectively. Thus, for each cooling period, the inequality (GST min,dk <GST min,d1 0) must be satisfied. N is the number of all the cooling periods in the vegetation phase. The higher the TDR value, the greater the temperature drop rate and the more severe the frost. Considering the relatively high temperature during the summer maize growing season and rare cold stress, the TDR was eliminated.
The standardized precipitation index (SPI) (Mckee et al 1993) is widely used to identify meteorological droughts and quantify the intensity and duration of droughts. The SPI has multiple time scales from 1-month to 24-months to depict different drought types. The SPI at a short time scale, such as 1 month, reflects short-term soil moisture conditions, while the SPI at time scale longer than 6 months monitors medium-and long-term precipitation trends, which relate to hydrologic regimes and reduce the covariation with vegetation vigor (Ji and Peters 2003). Previous studies demonstrated that 3-months SPI (SPI3) captures seasonal precipitation and agricultural drought risk well (Xu et al 2015). The response of crops to precipitation is not instantaneous but has a certain time lag. The NDVI and production anomaly of crops are most correlated with the seasonal cumulative moisture up to 3 months (Ji and Peters 2003, Patel et al 2007. Therefore, daily precipitation was integrated into monthly value to calculate the SPI3 by the method of Mckee et al (1993). Taking the moderate drought as a baseline (Mckee et al 1993, Zhang et al 2006, Zhai et al 2010, Lu et al 2020, we determined a meteorological drought event (DSPI3) when SPI3 was below −1. Similarly, an extreme wet event (WSPI3) occurred if SPI3 exceeded 1 (Zhai et al 2010, Xu et al 2015. The precipitation extremes (DSPI3, WSPI3) of wheat were measured by the SPI3 in May which represented the precipitation conditions from March to May covering the main wheat growing season. The SPI3 in September, which represented the precipitation conditions from July to September, was used to calculate precipitation extremes for maize.
In addition, the meteorological factors reflecting the mean climate conditions were also selected, including average growing season temperature (T mean ,°C), total precipitation during the growing season (P, mm) and in the previous season (preP, mm), average growing season solar radiation (SRD, MJ/(m 2 ·d)) and average growing season vapor pressure deficit (VPD, kpa) (Lobell and Burke 2010a, Lobell et al 2011b, Tao et al 2015, Tao et al 2017, Xiao et al 2018, Li et al 2019a. Among them, the pre-season accumulated precipitation for winter wheat included two variables: the total rainfall in summer maize season before wheat sowing (preP1, calculated from SOS to EOS of summer maize in the previous year) and the total rainfall from wheat sowing to green-up (preP2, calculated from a week after EOS of maize in the previous year to SOS of wheat). The preP for summer maize is the total rainfall in the wheat season before maize sowing (from SOS to EOS of wheat). Solar radiation (SRD) in each grid was estimated by sunshine hours and the Ångström-prescott equation (Zotarelli et al 2010). The average VPD during the growing season (SOS to EOS) was calculated by equation (4) (Zotarelli et al 2010, Lobell et al 2014 with mean temperature (T mean ) and minimum relative humidity (RH min ). Using minimum relative humidity rather than mean relative humidity was used to depict more severe water stress. We chose the averaged NDVI value during growing season (NDVI mean ) as the explained object because of its good approximation of crop growth status (Boken and Shaykewich 2002, Esquerdo et al 2011, Huang et al 2013, given that the yield data from statistics may not be reliable . In each grid, the daily NDVI data during the extracted wheat/maize growing season (SOS to EOS) were averaged to obtain the NDVI mean . The setting of the model considers the following three aspects. Firstly, the dependent variable was the NDVI mean of winter wheat or summer maize, and the explanatory variables contained climate and non-climate factors. Among them, the unmeasurable non-climate factors, such as the crop variety improvement, economic growth, soil type, etc, were controlled by two intercept terms. Secondly, irrigation and precipitation jointly constitute the water supply for crop growth, especially for winter wheat in the NCP. However, due to the unavailability of irrigation, the irrigation variable could not be used separately but was regarded as random error. Thus, the precipitation quantified variables (i.e., P, prep, DSPI3 and WSPI3) retained in the model may not represent the actual moisture condition, resulting in their attenuated impacts on the crop NDVI than under pure rain-fed cropland. Thirdly, necessary analyses were implemented. The unit root test was applied to check the stationarity of data (Jeločnik et al 2019), showing that all variables of winter wheat and summer maize were stationary. The Hausman test to determine the model types suggested that the fixed effects model was more appropriate (Jeločnik et al 2019). We calculated the correlation coefficients (tables S1 and S2) and the variance inflation factor (VIF) (table S3) of the variables to measure the multicollinearity among the variables. The correlation coefficients were low, and the VIF value of each explanatory variable also did not exceed 10, which was within a reasonable range. Neither analyses showed serious multicollinearity among the variables in the model (Xiao et al 2018, Jeločnik et al 2019, indicating that the influence of excessive collinearity on the final results could be ignored. Thus, the individual and time fixed effects panel regression model was determined to be the most reasonable equation.

Panel regression model
Equation (5) is for winter wheat.
, , In this model, NDVI mean i,t is the NDVI variable in grid i and year t. ln is the natural logarithm to reduce the heteroskedasticity of the dependent variable and make climate factors proportionally affect NDVI. α i is an intercept term with grid-fixed effects, representing the unobserved variables that are time-invariant but differ by grids (e.g., soil type, topography, and socio-economic characteristics). γ t is an intercept term with year-fixed effects, accounting for the unobserved variables which are time-varying but common in all grids (e.g., variety and technology improvement, national policy changes). β is the regression coefficient that represents the sensitivity of NDVI to corresponding climate variables. It represents the percentage variation of a dependent variable (β × 100(%)) caused by one unit increase of the corresponding climate variable keeping other explanatory variables constant. ε i,t is the random error term.
For summer maize, the panel regression model was similar to equation (5) without the term of TDR and only one variable of the pre-season total precipitation (preP) remained. The regression models of winter wheat and summer maize were fitted in the NCP.
In this study, we calculated the actual contribution of a climate variable to NDVI variation (

Assessment of phenology identification
The comparison between the phenological extraction results and the observed phenology at the agrometeorological stations is illustrated in figure 3. SOS, PEAK, and EOS for winter wheat corresponded its greenup, heading, and maturity dates, respectively. The R-square (R 2 ) was 0.96 (p<0.05), and the root mean square error (RMSE) was 10.3 days. For summer maize, SOS, PEAK, and EOS represented the seventh true leaf, heading, and maturity dates, respectively, with R 2 of 0.93 (p<0.05) and RMSE of 10.9 days. Variations in the three phenological dates and growing season length for winter wheat and summer maize are analyzed in the supplementary material (Text S2, figures S2, and S3).

Spatio-temporal changes of climate variables and the NDVI mean
We calculated the average values of climate variables and NDVI mean during winter wheat and summer maize growing seasons of all planting area grids in the NCP each year. The spatial patterns of the trends of climate variables and NDVI mean in each grid were also achieved. In the winter wheat season (figure 4), TDR significantly decreased with a rate of −0.01 (°C/d) year −1 , indicating an obvious weakening of spring frost over the past 30 years, while none of the trends of EDD, DSPI3 and WSPI3 reached a significance level. For mean climate factors, the growing season mean temperature showed a significant increasing trend at a rate of 0.03°C year −1 and the significant rise rate of VPD was 4.3e-03 kpa year −1 . No significant changes in solar radiation, total precipitation during growing season and in previous stages occurred. Spatially (figure 5), though the heat stress increased in most planting areas (accounting for 87.7%), the trend reached a significant level only in the extremely small areas scattered in the western and northern corners of the NCP (taking up 8.1%). The TDR concentrated in the central part of the NCP (only 2.2% of the planting area) rose significantly, while it declined significantly with an average rate of −0.04 (°C/d) year −1 in the areas evenly distributed in the northeast and south of the NCP (with a proportion of 19.9%). The tendencies of DSPI3 and WSPI3 in the NCP were not significant. Regarding the mean climate factors, the significantly increasing T mean was evenly distributed in most areas except for the northeastern part of the NCP (50.8% of the planting area), with an average increasing trend of 0.05°C year −1 . The SRD focused on the central and southern parts of the NCP decreasing significantly, with an average rate of −0.03 (MJ/(m 2 ·d)) year −1 . The areas where VPD increased significantly with an average rate of 5.8e-03 kpa year −1 were scattered, mainly in the northwest, central and south parts of the NCP. Areas with significant changes in P, preP1 and preP2 were rare and could be ignored.
In the summer maize season (figure 6), the heat stress increased at 4.96°C·h year −1 (p = 0.06). About extreme precipitation events, neither DSPI3 nor WSPI3 changed significantly. Except for the significant reduction of solar radiation at a rate of −0.03 (MJ/(m 2 ·d)) year −1 , the remaining mean climate variables exhibited insignificant increasing tendencies. In space (and figure 7), for climate extremes, 99.3% of the planting area experienced rising EDD, only 13.5% of the planting area dispersed in the southern NCP had a significant trend. The trends of DSPI3 and WSPI3 failed to reach a significant level in the entire NCP. For mean climate factors, the significantly increasing tendency of T mean was found mainly in the southeastern NCP. The SRD reduced significantly at an average rate of −0.05 (MJ/(m 2 ·d)) year −1 , which was evenly distributed in most areas except for the southeast. Only 6.0% of the planting area containing the west and the northeast corner of the NCP had a significantly increasing VPD. And the trends of P and preP were insignificant in the entire region.
In addition, the NDVI mean of winter wheat in the NCP significantly increased at a rate of 4.6e-03 year −1 from 1982 to 2015 (figure 4(k)). A significant increasing trend of NDVI mean was found throughout the study area ( figure 5(k)). Annual NDVI mean for summer maize showed an insignificant increasing trend (figure 6(i)). The spatial pattern of the significant growth of NDVI mean was relatively scattered, concentrated in the southern and northeastern regions (figure 7(i)).

The sensitivity of NDVI mean to climate variables
The sensitivities of the NDVI mean of winter wheat and summer maize to each climate variable reflected by the model coefficients (β) are shown in table 1. For winter wheat, though the coefficients of EDD and TDR showed that heat stress and spring frost had a positive relationship with wheat growth, the p-values of the T-statistic were too high to achieve statistical significance. The negative coefficient of meteorological droughts also did not reach  during 1982during -2015, the points corresponding to the years without meteorological drought (i.e., DSPI3 = −1) are not displayed. In subgraph (j), the points corresponding to the years without extremely wet event (i.e., WSPI3 = 1)) are not displayed. a significant level. Winter wheat growth responded negatively to extreme wet events, as one unit increase in WSPI3 significantly reduced the NDVI mean by 3.2% with a p-value close to 0.05 significance level. Besides, the NDVI mean showed higher sensitivity to extreme precipitation events than to extreme temperature events inferred by the higher absolute values of the coefficient of DSPI3/WSPI3 than that of EDD/TDR. Among mean climate factors, T mean and SRD significantly benefited NDVI mean , while VPD had the greatest negative impact on NDVI mean with the highest absolute coefficient value. The absolute values of their coefficients were noticeably greater than that of climate extremes. The total precipitation during the growing season (P) and in the previous seasons (preP1 and preP2) failed to show a significant impact on NDVI mean .
For summer maize, the coefficients of heat stress and meteorological droughts failed to reach statistical significance, and only extreme wet events caused significantly adverse impacts on NDVI mean for one unit increase in WSPI3 value would reduce NDVI mean by 2.0%. The impacts of mean climate variables on summer Figure 5. Spatial patterns of the trends of (a) mean temperature (T mean ), (b) total precipitation (P), (c) total precipitation in summer maize season before wheat sowing (preP1), (d) total precipitation from wheat sowing to green-up (preP2), (e) solar radiation (SRD), (f) vapor pressure deficit (VPD), (g) extreme degree days (EDD), (h) spring frost temperature-drop-rate (TDR), (i) meteorological droughts (DSPI3), (j) extreme wet events (WSPI3), and (k) NDVI mean during winter wheat growing season from 1982 to 2015. The statistical significance was p-value  0.05. maize growth were still greater than that of climate extremes. T mean and SRD significantly favored the NDVI mean , whereas VPD did not play a remarkably negative role on NDVI mean . The NDVI mean was also insensitive to P and preP.

The contributions of climate variables to NDVI mean during 1982-2015
From 1982 to 2015, we calculated the variation in the NDVI mean of winter wheat for all pixels, showing a median increase of 29.6% of NDVI mean . In figure 8(a), based on the median estimates, WSPI3 contributed the most to the change of NDVI mean over the past 34 years with a proportion of 0.3% among the climate extremes, and DSPI3 contributed the least (0.1%). The contributions of EDD and TDR were 0.2% and −0.1%, respectively. For mean climate variables, VPD, T mean , and SRD were predominant variables, leading to −6.6%, 5.2%, and −0.8% of the NDVI mean variation, respectively. The much higher absolute value of the contribution of the three mean climate variables also demonstrated that the variation of winter wheat growth was more restricted by mean climate conditions than by climate extremes. For summer maize, the NDVI mean increased by 3.6% during 1982-2015 based on the median estimate. All three climate extremes contributed negatively to the change of NDVI mean , and WSPI3 which played a significantly negative role contributed −0.1% of the NDVI mean change in the last 34 years. Moreover, T mean and SRD occupied a significant position for the variation. The maximum contribution to the NDVI mean variation was from T mean (a median value of 1.4%), followed by SRD (a median value of −1.2%) ( figure 8(b)).
We compared the contribution of the climate variables for each pixel and regarded the variable with the largest absolute contribution as the dominant climate factor. The result ( figure 9(a)) showed the leading role of the mean climate factors (i.e., T mean , VPD and SRD) in winter wheat growth. T mean was the dominant factor in the south region, while the areas where VPD made a leading contribution were dispersedly distributed in the northeast, northwest, central parts and the south margin of the NCP. SRD contributed most in a small area, focusing on the western corner of the NCP. For summer maize, T mean predominated in the southern and central parts of the NCP. The areas with the highest SRD contribution were located in the northern NCP. And EDD contributed the most in the planting areas scattered in the southern edge and northeast corner of the NCP ( figure 9(b)).

Discussion
5.1. The impacts of climate variables on the growth of winter wheat and summer maize Based on historical climatic data and remote sensing data, we revealed the changes of various climate factors in the NCP over the past 34 years and explored their impacts on the growth of winter wheat and summer maize using panel regression. Figure 6. Time series of averaged (a) mean temperature (T mean ), (b) total precipitation (P), (c) total precipitation in winter wheat season before maize sowing (preP), (d) solar radiation (SRD), (e) vapor pressure deficit (VPD), (f) extreme degree days (EDD), (g) meteorological droughts (DSPI3), (h) extreme wet events (WSPI3), and (i) NDVI mean during summer maize growing season of the planting area grids in the NCP during 1982-2015. Specially, in subgraph (g), the points corresponding to the years without meteorological drought (i.e., DSPI3 = −1) are not displayed. In subgraph (h), the points corresponding to the years without extremely wet event (i.e., WSPI3 = 1)) are not displayed.
Firstly, though previous studies concluded the adverse role of heat stress, we found that heat stress during the growing season did not significantly affect either the winter wheat or summer maize. Heat stress affects crops differently at various phenological stages. The pre-heading heat stress can decrease leaf elongation and stem Table 1. Model coefficients and determination coefficients (adjusted R 2 ) for the winter wheat and summer maize NDVI mean in response to climate variables in the panel regression model.

Model coefficient (β)
Winter wheat Summer maize Mean temperature (T mean ) 5.0E-02 (0.00) 2.3E-02 (0.00) Total precipitation (P) 1.0E-04 (0.35) −1.7E-05 (0.37) Total precipitation in maize season before wheat sowing (PreP1), or total precipitation in wheat season before maize sowing (  , we did not detect a significant impact of spring frost on winter wheat growth. A statistical analysis of reported extreme weather disasters and national cereal production loss data across the globe also indicated no significant production effects from extreme cold events (Lesk et al 2016). The possible reason is that spring frost damages young spikes, causes infertile spikelets and reduces the number of grains per ear (Wu et al 2014, Liu et al 2020b, resulting in yield loss, while it affects leaves and external morphology of wheat plant weakly, causing NDVI mean responds insignificantly. Thirdly, the meteorological droughts monitored by the three-month SPI failed to show a significant impact on either winter wheat or summer maize, which agreed with previous studies and could be due to the alleviation effect of widespread irrigation during the growing season (Wang et al 2018, Yu et al 2018a, 2018b, Lu et al 2020. Fourthly, the results showed that the extreme wet events harmed crop growth, especially for summer maize. The detrimental impact was also found in MäKinen et al (2018) that 76% of the wheat in Europe suffered yield penalties after experiencing heaving rainfall, and in Li et al (2019b) that the extreme wet conditions reduced maize yield on average by −16.6% relative to the yield trend in the United States. Considering that the SPI is a relative index, some years with high SPI values for wheat season and maize season were chosen, and daily precipitation in these seasons exceeded 40 mm, which was high enough to affect crop growth through physically damaging plants, destroying pollen, and excessive soil moisture restricting root  respiration (Li et al 2019b). In addition, wheat and maize cultivars have been renewed frequently over the past 30 years in China. For example, maize varieties with lodging resistance and high seed volume-weight were promoted extensively (Sun et al 2014). Cultivar renewal could improve crop resistance to extreme weather (Ci et al 2011). Compared with climate extremes, mean climate conditions, especially the average growing season temperature, vapor pressure deficit and solar radiation, have much greater impacts and contributions than climate extremes, indicating their dominance among climate factors for winter wheat and summer maize in the NCP. Our results are consistent with some previous studies. For example, Tao et al (2015) showed that both the increases in growing degree days above 0°C (GDD) and the decreases in solar radiation contributed more than high temperature (>34°C) degree days (HDD). Liu et al (2014) also demonstrated that wheat yield was more sensitive to average temperature than heat stress from 2000 to 2009 within our study area. Among the mean climate factors, we found that average temperature had a high positive impact and greatly contributed to winter wheat and summer maize growth. However, previous studies on the sensitivity of regional crop yields to temperature change failed to reach a consensus. For example, the correlation analysis of T mean and spring crop NDVI in the Huaihe plain showed their significant negative relationship with a coefficient of −0.508 (Li et al 2019a), seeing that increasing temperature can advance phenology, shorten the growing period and induce more pests and diseases (Piao et al 2010, Lobell et al 2011b, Chen et al 2018. While some studies found crop yield benefited from the warming trend. For example, Zhang and Yao (2013) analyzed that the average temperature benefited northcentral region for wheat yield and northeast region for maize yield in China counties during 1980-2008, considering the accelerating leaf photosynthesis and development rate under warming (Prasad et al 2008, Li et al 2019a. The inconsistency may be related to the temperature level, since crops respond nonlinearly to temperature Roberts 2009, Lobell et al 2011a), which means that warming favors crops until reaching its optimum value (Tao et al 2016, Zhao et al 2016, Li et al 2019a. The respective maximum T mean of winter wheat and summer maize during 1982-2015 in our study were 15.8°C and 26.9°C, respectively, showing they were still within the optimum temperature range (around 22°C for wheat (Porter and Gawith 1999), and 29°C for maize (Schlenker and Roberts 2009)) and may benefit wheat and maize. Besides, the growing season length was not shortened on a large scale, but prolonged, especially for winter wheat ( figure S2(d)), which may lead to the harmless effect of rising average temperature. The above-mentioned inconsistencies were also pointed out in Zhao et al (2016). The sensitivity of wheat yield to temperature was negative from crop models (yield decreased by 4.6% on average per°C) and warming field experiments (yield decreased by 2.8% on average per°C), but the statistical models showed a positive response as yield increased by 1.0% on average per°C in the NCP. Secondly, the highest absolute value of the vapor pressure deficit coefficient verified its considerable limitation on crops, especially winter wheat. The significantly negative impact of VPD agreed with previous studies (Lobell et al 2014, Hsiao et al 2019, Yuan et al 2019, considering that rising VPD triggered stomatal closure, leading to lower photosynthetic rates, consumption of carbohydrate reserves and impaired leaf area development. Yuan et al (2019) suggested that air temperature and VPD were the most important contributions to NDVI variability by machine learning model, and the analysis based on a coupled crop and soil model also demonstrated that increased VPD had a greater negative impact on yield than temperature (Hsiao et al 2019). Thirdly, the reducing trend of solar radiation, especially during the summer maize season, agreed with previous studies (Yang et al 2011, Qi et al 2015 and was mainly due to the increasing anthropogenic aerosols caused by industrialization (Wang and Yang 2014). The positive response of winter wheat and summer maize to solar radiation was consistent with previous studies (Tao et al 2016, Tao et al 2017. Though unfavorable effects of solar radiation were found in some studies (Tao et al 2012, Yin et al 2016, they attributed it to water stress caused by the coordination of solar radiation and high temperature. Good water supply in our study area could eliminate such negative effects. Moreover, the model coefficients of total precipitation during growing season and in the previous stage were insignificant, which was also found in previous studies Tao 2014, Xiao et al 2018). The outcome may be because extra irrigation supplemented precipitation and weakened the response to rainfall. Finally, although the NDVI mean of winter wheat and summer maize increased, especially winter wheat, how crop growth changes in the future is unknown and needs further research. For example, in the Guanzhong Plain of China, the future climate conditions with increasing temperature and precipitation benefit irrigated winter wheat yield but reduce summer maize yields. Appropriate management options, such as irrigation and adjusting planting dates, are required to combat climate change (Saddique et al 2020a(Saddique et al , 2020b.

Implication of the changes in NDVI mean for crop yield changes
Previous studies have pointed out a good correlation between NDVI and yield seeing that the photosynthetic capacity and growth status of crops measured by spectral-vegetation indices are directly related to ultimately yield and widely regarded the NDVI variable as a reasonable production predictor (Huang et al 2013). To evaluate how NDVI reflect the historical yield in our study area, we collected several city-level statistic winter wheat and summer maize yield data. We averaged the statistic yields and the NDVI mean of wheat/maize of all the planting area grids in the cities to obtain their regional averages.
As illustrated in figure 10, the NDVI mean of winter wheat agreed well with the statistical data, especially the clear upward trend of both the datasets since 1997. The correlation coefficient between the NDVI mean and wheat yield was 0.702 (p<0.05). The summer maize yield fluctuated obviously, with an increase by 56.3% from 1984 to 2015. The correlation coefficient for summer maize was 0.368 (p<0.05), and the agreement of the NDVI mean and maize yield was not as ideal as that of wheat. Li et al (2019a) compared the changes in maximum NDVI of the mature phase and the corresponding provincial statistical productions for spring and autumn crops in the Huaihe Plain which was similar to the NCP. They concluded that the correlation coefficient between NDVI and spring crop production was 0.857, and the correlation coefficient for autumn crops was lower (0.491), which was consistent with ours. Both the results showed that some uncertainties still existed between NDVI and statistical data, especially for summer maize. On the one hand, the uncertainty may stem from unreliable statistical data, which was highlighted in Liu et al (2020a) that the production data have potential misreporting problems affected by agricultural policies, insufficient supervision, and traditional statistical measure methods. On the other hand, the final yield is related to both the crop growth status and aboveground biomass monitored by NDVI and the harvest index (HI) which is defined as the ratio of grain yield to aboveground biomass. In the NCP, the average HI of summer maize is around 0.5 (Wang et al 2020) (lower than the HI of winter wheat: 0.62, ), which may lead to the poor accordance of NDVI mean and yield for maize. Therefore, though the NDVI mean is a good indicator for growth status and aboveground biomass, it may not be able to fully replace the changes in yield, especially summer maize. In view of the excellent performance of winter wheat, the impact of climate factors on NDVI revealed by our study provides a good guidance for its field management. While for summer maize, our conclusion may need to be more cautiously considered.

Limitations
Several limitations exist in this study and can be explored in future work. Firstly, a variety of extreme climate events and mean climate factors were considered in the panel regression model. Although the correlation matrix and variance inflation factor (VIF) test of the variables in our study did not indicate high multicollinearity and the combined use of spatial and time-series data also reduced the error caused by correlation, the effects of some climate variables were not completely independent and could be enhanced or weakened if combined with other variables (Jeločnik et al 2019, Li et al 2019b. For example, rising temperature without adequate water supply could be accompanied by an increase in VPD which can exert a comparable effect with temperature (Hsiao et al 2019). Besides, the four individual climate extremes and their impacts on crop growth were analyzed in the study, while the impacts of compound climate extremes were not taken into account, which has attracted more attention recently. For example, the concurrence of drought and heat induce insect outbreaks and damage natural systems (Zscheischler andSeneviratne 2017, Sun et al 2019). Compound soil moisture and heat events were more damaging than dry heat for maize in the USA (Haqiqi et al 2021). The severity of compound climate extremes can be considered in future research.

Conclusion
In the study, combining the long-series meteorological data and remote-sensing NDVI data, the impacts of climate changes on the growth of winter wheat and summer maize in the NCP from 1982 to 2015 were investigated by using the panel regression model. During the winter wheat season, spring frost decreased significantly, while the changes in other three climate extremes were insignificant changes. The growth of winter wheat over the past 34 years was not significantly sensitive to heat stress, spring frost, meteorological droughts, or extreme wet events, which may benefit from the compensatory effect of cultivar improvement, agricultural management practices and irrigation. Compared to climate extremes, wheat growth was much more sensitive to seasonal mean air temperature, solar radiation and vapor pressure deficit. The rising average temperature could accelerate the development rate, while the significantly increased vapor pressure deficit and insignificantly reduced solar radiation adversely affected leaf photosynthesis, contributing 5.2%, −6.6%, and −0.8% to the trend of wheat growth, respectively, which were much higher than that of climate extremes.
During the summer maize season, the extreme wet event was the only extreme climate factor exerting a significant impact on maize growth during 1982-2015, which could damage summer maize by restricting root respiration and destroying pollen. Heat stresses, which occurred frequently during the maize season, showed no significant impact. Meteorological droughts also played an insignificant role, indicating the complementarity of irrigation with insufficient precipitation. Among all climate factors, only solar radiation decreased significantly and summer maize responded significantly to average temperature and solar radiation. From 1982 to 2015, the extreme wet events contributed −0.1% of the changes in maize growth, which was less than the contributions of average air temperature (1.4%) and solar radiation (−1.2%). Our results reveal the important impacts of mean climate conditions on the trend of crop growth in dual-cropping systems and the adverse role of extreme wet events on summer maize in the NCP.