• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    Potential global distribution of the guava root-knot nematodeMeloidogyne enterolobii under different climate change scenarios using MaxEnt ecological niche modeling

    2023-07-17 09:42:46PANSongPENGDeliangLIYingmeiCHENZhijieZHAIYingyanLIUChenHONGBo
    Journal of Integrative Agriculture 2023年7期

    PAN Song,PENG De-liang,,LI Ying-mei,CHEN Zhi-jie,ZHAI Ying-yan,LIU Chen,HONG Bo#

    1 Shaanxi Key Laboratory of Plant Nematology,Bio-Agriculture Institute of Shaanxi,Shaanxi Academy of Sciences,Xi’an 710043,P.R.China

    2 State Key Laboratory for Biology of Plant Diseases and Insect Pests,Institute of Plant Protection,Chinese Academy of Agricultural Sciences,Beijing 100193,P.R.China

    Abstract In recent years, Meloidogyne enterolobii has emerged as a major parasitic nematode infesting many plants in tropical or subtropical areas. However,the regions of potential distribution and the main contributing environmental variables for this nematode are unclear. Under the current climate scenario,we predicted the potential geographic distributions of M. enterolobii worldwide and in China using a Maximum Entropy (MaxEnt) model with the occurrence data of this species. Furthermore,the potential distributions of M. enterolobii were projected under three future climate scenarios(BCC-CSM2-MR,CanESM5 and CNRM-CM6-1) for the periods 2050s and 2090s. Changes in the potential distribution were also predicted under different climate conditions. The results showed that highly suitable regions for M. enterolobii were concentrated in Africa,South America,Asia,and North America between latitudes 30°S to 30°N. Bio16 (precipitation of the wettest quarter),bio10 (mean temperature of the warmest quarter),and bio11 (mean temperature of the coldest quarter) were the variables contributing most in predicting potential distributions of M. enterolobii. In addition,the potential suitable areas for M. enterolobii will shift toward higher latitudes under future climate scenarios. This study provides a theoretical basis for controlling and managing this nematode.

    Keywords: Meloidogyne enterolobii,species distribution model,MaxEnt,climate change,future climate scenarios,centroid change

    1.lntroduction

    Root-knot nematodes (Meloidogynespp.),known for their wide geographical distribution and wide host range,have become one of the most threatening pathogenic nematodes in the world. One species ofMeloidogynespp.,Meloidogyne enterolobii,has emerged as a major parasitic nematode infesting many plants including various ornamentals,fruits trees,vegetable crops,wild plants,and weeds (Limaet al.2005;Carneiroet al.2006;Britoet al.2010;Poornimaet al.2016;Schwarzet al.2020).Meloidogyne enterolobiiis also known as the most aggressiveMeloidogynespecies because of its high reproduction rate and ability to overcome the resistance of some important crops,such as potato with theMhgene,tomato with theMi-1gene,bell pepper with theNgene,sweet pepper with theTabascogene,and soybean with theMir1gene (Berthouet al.2003;Britoet al.2007;Cetintaset al.2008). In a microplot experiment,tomato yield losses of up to 65% caused byM.enterolobiialone have been observed (Cetintaset al.2007). Many farmers might not even realize their fields have been infested withM.enterolobiiuntil the end of the season when crops are harvested and obvious galled root systems can be observed (Schwarz 2019). In heavily infested areas,cultivation may become unviable,such as guava cultivation in Brazil and India (Carneiroet al.2007;Singh 2020).

    Meloidogyne enterolobiiwas originally described on Pacara earpod trees (Enterolobium contortisiliquum) in Hainan Province of China in 1983 and has been reported in many countries (Yang and Eisenback 1983). Until recently,it was generally considered that the distribution ofM.enterolobiiwas restricted to regions with typical tropical or subtropical climatic conditions,such as Africa,South and Central America,Southeast Asia,and South Asia (EPPO 2022). In the United States,M.enterolobiiwas first reported in Puerto Rico and Florida. Since then,it has spread and has been reported in North and South Carolina (Britoet al.2004;Rutteret al.2019;Schwarzet al.2020). Apart from Hainan,M.enterolobiihas also been reported in other southern regions of China,including Guangdong,Guangxi,Fujian,Hunan,Yunnan,and Taiwan (Niuet al.2012;Wanget al.2015;Xiaoet al.2018;Lianget al.2020;Zhanget al.2020). In Brazil,M.enterolobiihas been detected in 18 states and is the most important pathogen for guava (Carneiroet al.2001;Souzaet al.2006;Limaet al.2007;de Almeidaet al.2008;Silvaet al.2008;Grothet al.2017;Silvaet al.2017;Luquiniet al.2019;Galbieriet al.2020). Apart from tropical and subtropical regions,M.enterolobiihas also been detected in commercial greenhouses in Switzerland,which is located in a temperate region (Kiewnicket al.2009).

    Climate change can influence the distribution and abundance of invasive pathogens directly and indirectly.Many studies have analyzed the invasion and distribution of invasive insects influenced by climate change,but few studies have reported on plant nematodes in this field. Global warming provides a potential opportunity forM.enterolobiito spread from low-latitude areas to highlatitude areas. Since it is probable thatM.enterolobiican survive in warmer areas and greenhouses throughout regions where no disease occurs,the risk of its establishment is high. Therefore,modeling the impact of climate change on the distribution ofM.enterolobiican provide vital information for controlling and managing the spread of this nematode.

    Some studies have applied species distribution models(SDMs) to predict the potential impact of invasive species to facilitate early warning and planning for future impacts(Ekesiet al.2016;Wei et al.2020). In this study,the current and future potential distribution ofM.enterolobiiwere estimated based on occurrence data using MaxEnt Software. In this study,we aimed to predict the trends of distribution and spread ofM.enterolobii,thereby providing a theoretical basis for the prevention and control of this nematode.

    2.Materials and methods

    2.1.Species occurrence data

    Occurrence data ofM.enterolobiiwere derived from databases of the European and Mediterranean Plant Protection Organization (EPPO,https://gd.eppo.int/),Center for Agriculture and Bioscience International (CABI,https://www.cabi.org/) and Global Biodiversity Information Facility (GBIF,https://www.gbif.org/). In addition,we also performed a literature search was performed in the Web of Science database using “Meloidogyne enterolobii”as a keyword and the Chinese literature database of China National Knowledge Infrastructure (CNKI,https://www.cnki.net/). Geo-referenced coordinates for each distribution sites were either obtained from the literature or using Google Earth Pro v7.3.4 (https://earth.google.com/web/). To check and minimize spatial biases,ENMTools 1.4 (Warrenet al.2010) was used to trim duplicate occurrences such that only one distribution point was present in each grid cell with a spatial resolution of 2.5 arc-minutes (approximately 4.5 km). In total,199 points forM.enterolobiiwere used in this study (Fig.1;Appendix A).

    2.2.Environmental variables

    Fig. 1 Global occurrence data of Meloidogyne enterolobii (red dot) for the MaxEnt modeling.

    Fig. 2 Receiver operating characteristic (ROC) curve generated by the MaxEnt model. The plot represents the sensitivity(true positive rate) and the specificity (false positive rate) of the model. The area under the ROC curve (AUC) represents the entire area underneath the ROC curve (red);the 95%confidence intervals are indicated in blue.

    Fig. 3 Jackknife test of variable importance for the MaxEnt model of M.enterolobii distribution. Light sea green,without a variable;blue,with only a single variable;red,with all variables.

    Fig. 4 Response curves of three environmental variables with the largest contribution. A,bio10 (mean temperature of warmest quarter). B,bio11 (mean temperature of coldest quarter). C,bio16 (precipitation of wettest quarter). The red curve represents the average of ten replicates;blue margins represent the standard deviation(±SD).

    Fig. 6 Global potential distribution map for Meloidogyne enterolobii under the future climate scenario. A,ssp126 in the year 2050.B,ssp126 in the year 2090. C,ssp585 in the year 2050. D,ssp585 in the year 2090. Gray,unsuitable habitat;green,low habitat suitability;orange,medium habitat suitability;red,high habitat suitability.

    Fig. 7 Potential distribution map for Meloidogyne enterolobii in China under the future climate scenario. A,ssp126 in the year 2050. B,ssp126 in the year 2090. C,ssp585 in the year 2050. D,ssp585 in the year 2090. White,unsuitable habitat;green,low habitat suitability;orange,medium habitat suitability;red,high habitat suitability.

    Fig. 8 Changes in the centroid of the potential distribution of Meloidogyne enterolobii in China. The black star represents the current centroid;dots of different color represent future centroids based on two emission scenarios (ssp126 and ssp585) in 2050 and 2090. The red and blue broken lines indicate the shift of centroid under ssp126 and ssp585 scenarios,respectively.

    Twenty-three environmental variables,including 19 bioclimatic variables,1 topographical variable,and 3 soil variables,were collected as the candidate variables for modeling. The topographical variable (altitude) and bioclimatic variables under current climate conditions(bio1-bio19,recorded from 1970 to 2000) were downloaded from the Worldclim Dataset v.2.1 (https://www.worldclim.org/data/worldclim21.html). Three soil variables (soil pH,soil organic carbon and soil texture)were extracted from the Harmonized World Soil Database v.1.2 (HWSD,http://www.fao.org/soils-portal/soil-survey/soil-maps-and-databases/harmonized-world-soildatabase-v12/en/). To avoid multicollinearity between these selected variables,we also calculated Pearson’s correlation coefficients using SDMtoolbox 2.4 (http://www.sdmtoolbox.org/) in ArcGIS 10.2 (ESRI,Redlands,CA,USA) to remove highly correlated variables (r≥|0.8|)(Appendix B;Brownet al.2017). Finally,12 variables(including bio2,bio3,bio10,bio11,bio16,bio17,bio18,bio19,altitude,topsoil pH,topsoil organic carbon,and topsoil texture) were used to predict the distribution ofM.enterolobii.

    To predict the global potential distribution ofM.enterolobiiunder future climate conditions,8 selected bioclimatic variables available from three types of Global Climate Models (GCMs) including Beijing Climate Center Climate System Model version 2 (BCC-CSM2-MR),Canadian Earth System Model version 5 (CanESM5) and Centre National de Recherches Météorologiques Climate Model version 6 (CNRM-CM6-1) for 2050 (average for year 2041-2060) and 2090 (average for year 2081-2100),as well as 2 greenhouse gas emission scenarios including ssp126 (low greenhouse gas emission scenario)and ssp585 (high greenhouse gas emission scenario)according to the Shared Socioeconomic Pathways (SSP)under Phase 6 of the Coupled Model Intercomparison Project were obtained from the Worldclim Dataset. All of the current and future environmental variables were resampled to a 2.5-min spatial resolution using ArcGIS.To reduce the variations among different GCMs,a potential habitat suitability map forM.enterolobiiwas finally generated by averaging the projections of the three GCMs under all future climate scenarios. In addition,SDMtoolbox was used to analyze the centroid position changes of the total suitable habitat area ofM.enterolobiiin China under current and future climate scenarios.

    2.3.Optimization of model parameters

    A well-chosen set of model parameters can enhance the performance of the MaxEnt model. Recent studies have demonstrated that the default settings of MaxEnt may produce overfitted models and are not always appropriate for species distribution modeling (Merowet al.2013;Radosavijevic and Anderson 2014). Thus,to avoid overfitting and produce the optimal model forM.enterolobii,the kuenm package (Coboset al.2019) was used in R 4.1.0 (https://www.r-project.org/) to compare various combinations of the two most critical parameters,feature classes and regularization multiplier,to choose the best combination. MaxEnt contains 5 feature classes: linear(L),quadratic (Q),product (P),threshold (T),and hinge(H);thus,there were 31 combinations of feature classes in total. Moreover,8 values of the regularization multiplier varied from 0.5 to 4 with an interval of 0.5. Altogether,248 candidate models with all combinations of 8 regularization multiplier settings,31 feature class combinations,and a set of 12 environmental variables were evaluated.

    Model parameters were optimized using the functionkuenm_cevalbased on statistical significance (partial receiver operating characteristic,ROC,with 500 iterations),predictive ability (omission rate,OR),and model complexity (AICc) (Petersonet al.2008;Zhuet al.2017). According to the “OR_AICc” criterion,the significant candidate models below the threshold of omission rate (for example,≤0.05 when possible) and with the lowest ΔAICc values (≤2) were considered the final models with the optimal model parameters (Coboset al.2019).

    2.4.MaxEnt modeling

    We used MaxEnt v.3.4.1 (https://biodiversityinformatics.amnh.org/open_source/maxent/) in this study to model the potential suitable habitats forM.enterolobii. In this study,75% of occurrences were randomly selected for model training,and the remaining occurrences (25%) were used for model testing. The analysis included 10 replicates and a maximum number of 5 000 iterations involving the abovementioned 199 distribution points forM.enterolobiiand 12 environmental variables. The percent contribution,permutation importance,as well as jackknife test were used to determine each environmental variable’s relative contribution to the models (Phillipset al.2017). When the current distribution model forM.enterolobiiwas built,we projected it to 2 future periods (2050 and 2090) in MaxEnt.

    The logistic output was used to generate a map in raster format with species suitability values between 0 and 1. Fixed cumulative value 5 logistic threshold(FCV5) was used to define suitable and unsuitable habitats forM.enterolobii,and the areas with an output value below FCV5 were regarded as unsuitable habitats(Weiet al.2018). In ArcGIS,the Natural Breaks (Jenks)method was used to reclassify the suitable habitats forM.enterolobiiinto 3 levels: low habitat suitability,medium habitat suitability,and high habitat suitability. The areas of suitable habitats at different levels were calculated and distribution maps were generated using ArcGIS.Moreover,to render every model comparable in space and time series and simplify the analysis process,future climate scenario maps were also prepared using the same criteria as for the current climate scenario.

    2.5.Model evaluation

    The area under the ROC curve (AUC) was used to estimate the performance of the model (Swets 1988).The AUC value ranged from 0 to 1;AUC<0.5 indicates a random prediction,0.5≤AUC<0.7 indicates poor model performance,0.7≤AUC≤0.9 indicates moderate performance,and AUC>0.9 indicates high performance(Araujoet al.2005;Petersonet al.2011).

    3.Results

    3.1.Model performance for M.enterolobii

    Among the 248 candidate models,all models were statistically significant and only one model met OR and AICc criteria (delta AICc=0) (Appendices C and D). Thus,the model M_1.5_F_lt_Set1 (regularization multiplier=1.5,feature class combinations=L and T) was chosen in the final MaxEnt settings.

    The MaxEnt model performance forM.enterolobiiwas better than random with an average test AUC value of 0.939±0.017 (Fig.2),indicating that the model had an excellent performance. The global habitat areas forM.enterolobiiwere divided into four levels: (1) <0.095,which indicated an unsuitable habitat;(2) 0.095-0.228,which indicated low habitat suitability;(3) 0.228-0.450,which indicated medium habitat suitability;and (4) >0.450,which indicated high habitat suitability.

    3.2.Variable importance

    The percent contribution and permutation importance of the 12 environmental variables are listed in Table 1. The results showed that bio16 (precipitation of the wettest quarter) was the variable contributing the most. Its contribution rate and permutation importance reached 39.6 and 39.1%,respectively,indicating that bio16 was the main rainfall factor affecting the distribution ofM.enterolobii.In addition,the percent contribution of bio11 (mean temperature of the coldest quarter) reached 17.9% and the permutation importance of bio10 (mean temperature of the warmest quarter) reached 19.7%,indicating that bio11 and bio10 were the main temperature factors.

    Table 1 Percentage contribution of 12 environmental variables

    The results of the jackknife test yielded the relative importance of each variable (Fig.3). Among these variables,bio16,bio10,and bio11 were the leading three variables in terms of regularized training gain. Thus,the main variables for predicting the potential distribution ofM.enterolobiiwere bio10,bio11,and bio16. The remaining environmental variables had relatively low influence.

    According to response curves for the three environmental variables contributing most to the prediction ofM.enterolobiidistribution (Fig.4),the average ranges of each variable associated with high habitat suitability(probability>0.450) were 25.4-29.0°C for bio10,12.3-22.8°C for bio11,and 360-1 170 mm for bio16.

    3.3.Potential distribution of M.enterolobii under current climate scenario

    The global potential distribution ofM.enterolobiibased on current environmental variables is shown in Fig.5-A. The results showed that the regions suitable forM.enterolobiiwere distributed on all continents except Antarctica,and most suitable areas were concentrated in Africa,South America,Asia,and North America between latitudes 30°S and 30°N. The global land area of potential suitable habitat forM.enterolobiiwas approximately 5 147.28×104km2,42.9% of which (2 209.69×104km2) comprised regions with low suitability,37.1% of which (1 909.02×104km2) had medium suitability,and 20.0% (1 028.57×104km2) were highly suitable regions (Table 2).

    Table 2 Area with suitability for Meloidogyne enterolobii under different climate scenarios

    In addition,the distribution pattern ofM.enterolobiiwas significantly different on different continents. Africa had the largest suitable area forM.enterolobii(Fig.5-B),and areas with medium and high habitat suitability were mainly distributed in western,central,and southern Africa.In South America,the areas suitable forM.enterolobiicovered all 12 countries. Highly suitable regions forM.enterolobiiwere distributed in Brazil,Argentina,Paraguay,Bolivia,Venezuela,and Guyana. In Asia,the areas suitable forM.enterolobiiwere primarily distributed in East Asia (China,Japan,and South Korea),South Asia(India,Nepal,Bangladesh,and Sri Lanka),and Southeast Asia (Burma,Thailand,Vietnam,Cambodia,and the Philippines). In North America,the areas of potential distribution mainly included the southeast region of the United States,Mexico,Honduras,Nicaragua,and Cuba. In Oceania,the areas suitable forM.enterolobiiwere mainly distributed in northern and eastern Australia. Europe had the fewest suitable areas,and the total estimated suitable habitat area forM.enterolobiiwas approximately 312.09×104km2(Fig.5-B). Except for a few small areas of medium habitat suitability located in France,Portugal,and Italy,the remaining regions were regarded as areas of low habitat suitability or unsuitable forM.enterolobii.

    Under the current climate scenario,the area of potential suitable habitat forM.enterolobiiin China was approximately 300.14×104km2,accounting for 31.26%of the national area (Table 2). Highly suitable regions forM.enterolobiicovered 13.34% (128.11×104km2) of the total area of China and were mostly distributed in southern regions with latitudes lower than 35°N,including Hainan,Guangxi,Guangdong,and Fujian,most parts of Yunnan,Guizhou,Chongqing,Hunan,Jiangxi,Taiwan,eastern Sichuan,as well as small parts of Tibet and Zhejiang.Areas with medium habitat suitability accounted for 13.98%of the national area,covering 134.19×104km2. Areas with low suitability forM.enterolobiicovered 10.51% of the total area of China,totaling 37.84×104km2(Fig.5-C).

    3.4.Potential distribution of M.enterolobii under future climate scenarios

    The global potential distribution ofM.enterolobiibased on two emission scenarios (ssp126 and ssp585) in 2050 and 2090 is shown in Fig.6 and Table 2. Compared with the current climate scenario,the area of suitable habitat forM.enterolobiiin future climate scenarios increased to varying degrees. Under the ssp126 scenario,it is estimated that by 2050,the total area of suitable habitat will be 5 484.35×104km2,which is increased by 6.55%compared with the current area of suitable habitat. The area of highly suitable habitat will expand to 1 156.09×104km2,which is increased by 12.40% (Fig.6-A). In contrast to the current distribution,by 2090,the total suitable area and highly suitable area will increase by 7.34 and 12.54%,respectively,which exhibited a slightly increasing trend compared with that of 2050 (Fig.6-B). Under the ssp585 scenario,it is estimated that by 2050,the total suitable area and highly suitable area ofM.enterolobiiwill expand to 5 586.98×104and 1 212.99×104km2,which is increased by 8.54 and 17.93%,respectively (Fig.6-C). In contrast to the current distribution,by 2090,the total suitable area will increase by 18.57%,and the highly suitable area will increase by 4.11%. In that period,the total area of suitable habitat will be 6 103.06×104km2,which represents a significant increase compared with that of 2050. However,the area of highly suitable habitat in 2090 was smaller than that in 2050 for the ssp585 scenario(Fig.6-D). In addition,the highly suitable regions forM.enterolobiitended to spread from low latitudes to high latitudes. North America except for the southeast region of the United States,North and West Europe,Saharan Africa,Middle East,Central Asia,Russia,Mongolia,northern China,and southern Australia remained regions that had an unsuitable habitat forM.enterolobii(Fig.6).

    The potential distribution ofM.enterolobiiin China under the ssp126 and ssp585 scenarios is shown in Fig.7 and Table 2. Compared with the current climate scenario,the total suitable area and the area of highly suitable habitat forM.enterolobiiincreased,but the area of moderately suitable habitat decreased under two future climate scenarios. Under the ssp126 scenario,in the 2050s,the total area of suitable habitat will be 329.92×104km2,which is increased by 9.92%. The area of highly suitable habitat will expand to 227.75×104km2,which is increased by 73.87%,and some provinces,such as Hubei,Anhui,the southern parts of Henan and Jiangsu,will become highly suitable regions forM.enterolobiicompared with the current area of suitable habitat (Fig.7-A). By 2090,the total suitable area and highly suitable area will increase by 11.93 and 76.88%,respectively (Fig.7-B). Under the ssp585 scenario,the area of suitable habitat forM.enterolobiifurther increased. In the 2050s,the total suitable area and highly suitable area ofM.enterolobiiwill expand to 344.03×104and 228.14×104km2,which is increased by 14.62 and 78.08%,respectively (Fig.7-C). By 2090,the total area of suitable habitat will be 399.43×104km2,which is increased by 33.08% compared with the current area of suitable habitat. In that period,the eastern and southern parts of Shandong will become highly suitable regions forM.enterolobii,whereas the area of highly suitable habitat in Anhui,Hubei,Hunan,and Jiangxi will decrease compared with that of 2050. In addition,Central Shaanxi,Southern Shanxi,Eastern Gansu,and Southwestern Hebei will become moderately suitable regions forM.enterolobii.However,Xinjiang,Qinghai,Inner Mongolia,Heilongjiang,and Jilin remained regions of unsuitable habitat forM.enterolobii(Fig.7-D).

    3.5.Centroid changes in potential distribution

    The centroids of the habitat ofM.enterolobiiin China under different climate scenarios are shown in Fig.8.For the range ofM.enterolobii,the current centroid was located in Guangxi (110.589°E,24.921°N). Under the ssp126 scenario,the centroid of the suitable area is expected to shift northeastward into Shaoyang City,Hunan Province (111.185°E,27.052°N) during the 2050s,and to 111.305°E,27.207°N during the 2090s. Under the ssp585 scenario,the centroid of the future suitable area will be located at a northeastern position (110.886°E,27.188°N) during the 2050s but shifts in the northwestern direction to Huaihua City (110.347°E,27.442°N) during the 2090s. In general,the centroid of the potential suitable area ofM.enterolobiiwill migrate to higher latitude areas under the future climate scenarios.

    4.Discussion

    Compared with other niche models,MaxEnt has several advantages,such as easy operation,low sample number requirement,low uncertainty,and high accuracy. Thus,MaxEnt has been widely used in ecology,conservation biology,and invasive biology,especially in the prediction of potential distribution of invasive alien species in recent years (Phillipset al.2006;Liet al.2020). An appropriate model complexity from MaxEnt is critical to obtain a better-performing model in species distribution prediction(Merowet al.2013). However,the default parameters of MaxEnt might lead to an increase in model complexity and yield over-fitted prediction models,which result in reduced cross-temporal transferability and predictive credibility of MaxEnt models for forecasting species distributions under future climate change (Warren and Seifert 2011;Moreno-Amatet al.2015). In this study,two main parameters,feature classes and regularization multiplier,known to influence the prediction results by MaxEnt,were optimized to constrain the complexity of the model.The results indicated that a regularization multiplier=1.5 and feature class combinations=LT were selected as the optimal parameter set,and a model with the lowest AICc was used for the potential distribution prediction ofM.enterolobiiunder different climate scenarios. The estimated test AUC value of 0.939 indicated an excellent predictive accuracy of the model.

    Recently,bioclimatic variables and topographic variables have emerged as powerful tools for studying biological invasion and predicting potential global and regional distributions and habitat suitability (Jianget al.2018;Jinget al.2020;Weiet al.2020). However,different from insects,the life cycle of root-knot nematode is completed in the soil. Therefore,soil factors might significantly affect the distribution ofM.enterolobii. In this study,23 environmental variables (19 bioclimatic variables,1 topographical variable,and 3 soil variables) were obtained and their correlation analyses were performed using Pearson’s correlation.Finally,12 variables were selected for the distribution prediction ofM.enterolobii. Based on the jackknife test,the permutation importance,and contribution rate of these variables calculated using the MaxEnt model,the predicted results indicated that the mean temperature of the warmest quarter (bio10),mean temperature of the coldest quarter(bio11),and precipitation of the wettest quarter (bio16)were the main environmental variables affecting the distribution and regional changes in suitable habitat forM.enterolobii. Therefore,its distribution is determined by both temperature and precipitation factors. By contrast,other environmental factors,such as altitude,topsoil pH,topsoil organic carbon,and topsoil texture,had less impact on the distribution ofM.enterolobii. According to previous studies,temperature and humidity factors had significant effects on the survival of nematodes,with different temperatures and humidity required for the development and reproduction of different root-knot nematode species(Charwatet al.2002;Dávila-Negrón and Dickson 2013;Wuet al.2018). Compared with other species,M.enterolobiiwere more suited for reproducing and surviving in a high-temperature environment,and exposure to lower temperature decreased root penetration and slowed the developmental cycle ofM.enterolobii(Vellosoet al.2022). Rainfall also had a significant effect on nematode abundance and diversity. A significant trend for increasedMeloidogynespp.diversity and greater prevalence with increasing rainfall has been observed (Jordaanet al.1989;Fleminget al.2016). This is consistent with the characteristics of high temperature and precipitation in highly suitable habitats forM.enterolobiifound in this study(Fig.4). When temperature and precipitation increased,the temperature and humidity of the soil also increased accordingly,which was more conducive to the survival ofM.enterolobii.Conversely,cold and arid environments wereunfavorable for the survival ofM.enterolobii.

    In this study,the global potential distribution ofM.enterolobiiwas first simulated using MaxEnt. The potential global distribution ofM.enterolobiiunder the current climate scenario was concentrated in Africa,South America,Asia,and North America between latitudes from 30°S to 30°N,which is consistent with the occurrence data ofM.enterolobii. The suitable areas forM.enterolobiiwere mostly tropical and subtropical regions with large quantities of precipitation and high temperatures,which not only are more conducive to the reproduction and survival ofM.enterolobiibut alsoallow for greater plant diversity,which is more beneficial for the distribution and spread ofM.enterolobii. On the contrary,although North Africa is a tropical region,it was not suitable for the distribution ofM.enterolobiiowing to little rainfall and an extremely dry climate. It should be noted that in European countries such as France,Belgium,Switzerland,and Portugal,M.enterolobiiis mainly distributed in greenhouse environments,which differs from local field environments(Kiewnicket al.2008). To a certain degree,this might affect the accuracy of model prediction. The area of potential suitable habitat forM.enterolobiiin China was located south of 35°N. The reason for this phenomenon might be that the field environment is too cold for the survival ofM.enterolobiiin the winter. In some high northern latitudes where the occurrence ofM.enterolobiihas been reported,such as France and Shaanxi (China),the nematodes could not be detected one or two years after first being detected,which indicated that these areas might not be suitable for the distribution ofM.enterolobii(Santoset al.2019).

    The potential distribution ofM.enterolobiiunder future climate scenarios indicated that from now to 2090s,the suitable area ofM.enterolobiishowed a trend of gradual expansion and movement toward higher latitudes,which might be related to the trend of global warming. These results are consistent with some other invasive pests,such asDysmicoccus neobrevipes(Weiet al.2020) andCulex pipiens pallens(Liuet al.2020). Based on previous studies,under the medium greenhouse gas emission scenario (RCP4.5),the average temperature will rise by 3.0-5.2°C,and annual precipitation will increase by 4.2-6.1% in northern China by the end of the 21st century(Moet al.2018). These results suggested that compared with the current period,a higher temperature and humidity environment in northern China in the future will provide an ideal habitat forM.enterolobii. Therefore,the potential suitable area ofM.enterolobiiwill shift toward the northern part of China under future climate scenarios.

    In this study,SDMToolbox Software was used to analyze changes in the centroid of the habitat ofM.enterolobiiin China under different climate scenarios. Currently,the centroid is located in Guilin City,in Guangxi Zhuang Autonomous Region (110.589°E,24.921°N). In the 2050s and 2090s,the centroid of the total suitable area ofM.enterolobiishowed a trend of migration to Hunan Province under two different greenhouse gas emission scenarios (ssp126 and ssp585),which validated the results that climate warming would lead to the spread ofM.enterolobiifrom low-latitude regions to high-latitude regions. To better predict changes in the distribution ofM.enterolobiiin the future,it is necessary to analyze the ecological and physiological mechanisms of its adaptation to environmental factors,such as temperature and humidity.

    Considering the wide global distribution of the host andthe high level of agricultural trade between countries,there is a high risk ofM.enterolobiitransmission to non-suitable areas through susceptible plant materials(e.g.,parts of the roots or tuber stems),although the global distribution area of these nematodes is mainly concentrated in tropical and subtropical regions at present. In 2010,M.enterolobiiwasadded to the European and Mediterranean Plant Protection Organization A2 Alert list and recommended for regulation as quarantine pests in some states of the United States(USDA and PCIT 2014;EPPO 2022). A highly specific and sensitive diagnostic method to detectM.enterolobiiis necessary for monitoring the transmission and preventing the spread of this highly destructive nematode. To control the spread of this nematode species,it is essential to cultivate non-infected seedlings and the infected plant material,infested soil,and contaminated farm-equipment from infested fields withM.enterolobiishould not be moved to non-infested regions (Schwarzet al.2020).Meanwhile,several management strategies,such as soil solarization,crop rotation with non-host crops,biological control,new natural resistance sources,and fumigant and non-fumigant nematicide have been employed to control the spread ofM.enterolobiiworldwide (Zasadaet al.2010;Carrillo-Fasioet al.2020).

    As an endoparasite,M.enterolobiifeeds and matures to the adult stage of its life cycle fully inside host plant tissue (Elling 2013;Sureshet al.2019). Therefore,biological factors cannot be ignored in predicting the distribution ofM.enterolobii.In this study,the model used considered only the impact caused by abiotic factors but not that caused by biological factors,which may lead to predicting wider suitable areas ofM.enterolobii.It will be necessary to consider the impact of biological factors,such as the distribution ofM.enterolobiihost and different crop planting modes in our subsequent research.

    5.Conclusion

    In the present study,maps of the potential distribution ofM.enterolobiiin China and worldwide based on current and future climate scenarios are provided for the first time. The range of highly suitable habitats of this nematode increased and shifted toward higher latitudes under future climate scenarios compared with the current climate scenario. Therefore,the northern part of China is considered to have a risk of invasion in the future. We sought to predict the trends of distribution and spread ofM.enterolobii,which provides a theoretical basis for the prevention and control of this species.

    Acknowledgements

    This work was supported by the Key R&D Project of Shaanxi Province,China (2020ZDLNY07-06) and the Science and Technology Program of Shaanxi Academy of Sciences (2022k-11). We thank LetPub (www.letpub.com) for its linguistic assistance during the preparation of this manuscript.

    Declaration of competing interest

    The authors declare that they have no conflict of interest.

    Appendicesassociated with this paper are available on https://doi.org/10.1016/j.jia.2023.06.022

    国产av码专区亚洲av| 日韩av在线免费看完整版不卡| 精品视频人人做人人爽| 亚洲丝袜综合中文字幕| 免费看不卡的av| 黄色配什么色好看| 亚洲欧美精品自产自拍| 日韩欧美 国产精品| 男女免费视频国产| 国产精品无大码| 欧美bdsm另类| 亚洲一区二区三区欧美精品| 国产成人精品久久久久久| 一级av片app| av女优亚洲男人天堂| 久久精品久久精品一区二区三区| 免费不卡的大黄色大毛片视频在线观看| 人妻夜夜爽99麻豆av| 亚洲中文av在线| 免费人成在线观看视频色| 久久精品国产亚洲av天美| 最近2019中文字幕mv第一页| 国产成人免费无遮挡视频| 黑人高潮一二区| 日产精品乱码卡一卡2卡三| 人体艺术视频欧美日本| 国产探花极品一区二区| 国产一区有黄有色的免费视频| 国产日韩欧美视频二区| 18禁在线无遮挡免费观看视频| 亚洲欧洲精品一区二区精品久久久 | 成人二区视频| 国语对白做爰xxxⅹ性视频网站| 青春草亚洲视频在线观看| 亚洲精品乱久久久久久| 日韩伦理黄色片| 国产永久视频网站| 欧美日韩综合久久久久久| 中文在线观看免费www的网站| 丝袜脚勾引网站| 国产成人免费无遮挡视频| 丰满乱子伦码专区| 午夜激情久久久久久久| 最新中文字幕久久久久| 精品人妻偷拍中文字幕| 久久久久国产网址| 中文字幕久久专区| 亚洲av综合色区一区| 一级毛片我不卡| 嫩草影院新地址| 九色成人免费人妻av| 菩萨蛮人人尽说江南好唐韦庄| 国产在线免费精品| 青春草视频在线免费观看| 黑人巨大精品欧美一区二区蜜桃 | 亚洲精品aⅴ在线观看| 亚洲av二区三区四区| 熟女av电影| 成人美女网站在线观看视频| 狠狠精品人妻久久久久久综合| 人妻少妇偷人精品九色| 在线观看免费高清a一片| 久久6这里有精品| 亚洲,一卡二卡三卡| 欧美日韩国产mv在线观看视频| 美女cb高潮喷水在线观看| 美女cb高潮喷水在线观看| 国产精品一区www在线观看| 好男人视频免费观看在线| av在线播放精品| 青春草亚洲视频在线观看| 狂野欧美激情性xxxx在线观看| 大片电影免费在线观看免费| xxx大片免费视频| 99久久中文字幕三级久久日本| 亚洲伊人久久精品综合| 九色成人免费人妻av| 国产免费视频播放在线视频| 少妇人妻精品综合一区二区| 丝袜喷水一区| 久久久久国产网址| 免费观看a级毛片全部| 亚洲内射少妇av| 免费观看的影片在线观看| 久久久久久久国产电影| 一个人看视频在线观看www免费| 国内少妇人妻偷人精品xxx网站| 亚洲熟女精品中文字幕| 777米奇影视久久| 夫妻午夜视频| 好男人视频免费观看在线| 精品一区二区三区视频在线| 亚洲精品乱码久久久v下载方式| 在线精品无人区一区二区三| 男人舔奶头视频| 成人特级av手机在线观看| 能在线免费看毛片的网站| 亚洲人与动物交配视频| 色视频www国产| 内地一区二区视频在线| 乱码一卡2卡4卡精品| 久久这里有精品视频免费| 成人免费观看视频高清| 久久久午夜欧美精品| 国内揄拍国产精品人妻在线| 夜夜看夜夜爽夜夜摸| freevideosex欧美| 久久久久久久久大av| 天堂8中文在线网| 麻豆成人午夜福利视频| 欧美成人午夜免费资源| 国产av一区二区精品久久| 亚洲精品乱码久久久久久按摩| 91久久精品国产一区二区三区| 青春草国产在线视频| 久久av网站| 久久精品国产自在天天线| 精品卡一卡二卡四卡免费| 亚洲av在线观看美女高潮| 99视频精品全部免费 在线| 国产欧美另类精品又又久久亚洲欧美| 夜夜骑夜夜射夜夜干| 中文字幕制服av| 精品一区二区三卡| 99九九在线精品视频 | 观看免费一级毛片| 午夜免费观看性视频| 亚洲精品国产成人久久av| 2021少妇久久久久久久久久久| 纯流量卡能插随身wifi吗| 99re6热这里在线精品视频| 中文天堂在线官网| 国产欧美日韩综合在线一区二区 | 一级毛片久久久久久久久女| 九色成人免费人妻av| 少妇 在线观看| 久久毛片免费看一区二区三区| 成人毛片60女人毛片免费| 亚洲国产精品专区欧美| 在线观看美女被高潮喷水网站| 亚洲国产色片| 精品亚洲成国产av| 久久亚洲国产成人精品v| 亚洲av二区三区四区| 日本av免费视频播放| 国产精品欧美亚洲77777| 一区二区三区乱码不卡18| 亚洲人成网站在线播| 免费看光身美女| 国产成人91sexporn| 成年人午夜在线观看视频| 久久久欧美国产精品| 如日韩欧美国产精品一区二区三区 | 26uuu在线亚洲综合色| 一级毛片我不卡| 人妻 亚洲 视频| 2022亚洲国产成人精品| 精品久久久精品久久久| 久热这里只有精品99| 日产精品乱码卡一卡2卡三| av国产精品久久久久影院| 在线观看人妻少妇| 国产黄色免费在线视频| 亚洲精品国产av成人精品| 亚洲丝袜综合中文字幕| 国产真实伦视频高清在线观看| av国产精品久久久久影院| 免费黄频网站在线观看国产| 亚洲精品日韩在线中文字幕| 卡戴珊不雅视频在线播放| 中文在线观看免费www的网站| 男的添女的下面高潮视频| 国产亚洲最大av| 如何舔出高潮| 男女边摸边吃奶| 综合色丁香网| 少妇熟女欧美另类| 中文字幕制服av| 国产黄片视频在线免费观看| 最后的刺客免费高清国语| 国产无遮挡羞羞视频在线观看| 赤兔流量卡办理| av福利片在线观看| 国产 一区精品| 七月丁香在线播放| 狂野欧美激情性xxxx在线观看| 精品一区二区三区视频在线| 伊人久久国产一区二区| 免费看光身美女| 亚洲美女视频黄频| 欧美日韩av久久| 99国产精品免费福利视频| 日本vs欧美在线观看视频 | 夜夜看夜夜爽夜夜摸| 男女国产视频网站| 亚洲图色成人| 精品人妻偷拍中文字幕| 日韩欧美精品免费久久| 欧美人与善性xxx| 桃花免费在线播放| 日本黄色片子视频| 中文天堂在线官网| 嫩草影院入口| 亚洲,欧美,日韩| 免费不卡的大黄色大毛片视频在线观看| 一级黄片播放器| 国产欧美亚洲国产| 成人特级av手机在线观看| 中文字幕亚洲精品专区| 91午夜精品亚洲一区二区三区| 国精品久久久久久国模美| 午夜影院在线不卡| 日本av免费视频播放| 性色avwww在线观看| 免费看不卡的av| 亚洲欧美成人综合另类久久久| 中文乱码字字幕精品一区二区三区| 日日摸夜夜添夜夜添av毛片| 自拍欧美九色日韩亚洲蝌蚪91 | 菩萨蛮人人尽说江南好唐韦庄| 熟女人妻精品中文字幕| av天堂中文字幕网| 久久99热这里只频精品6学生| 国产高清三级在线| 免费播放大片免费观看视频在线观看| 人妻系列 视频| 欧美区成人在线视频| 国产视频首页在线观看| 国产一区有黄有色的免费视频| av播播在线观看一区| 国语对白做爰xxxⅹ性视频网站| kizo精华| xxx大片免费视频| 成人美女网站在线观看视频| 日本av免费视频播放| 99久久精品国产国产毛片| 欧美+日韩+精品| 各种免费的搞黄视频| 韩国高清视频一区二区三区| 国产av精品麻豆| 2018国产大陆天天弄谢| 狂野欧美激情性bbbbbb| 青春草视频在线免费观看| 中文字幕av电影在线播放| 九九爱精品视频在线观看| 黄色日韩在线| 观看av在线不卡| 一级毛片久久久久久久久女| 一区在线观看完整版| 最新的欧美精品一区二区| 亚洲精品国产成人久久av| 日本欧美国产在线视频| 久久久久久久久久久丰满| 黄色一级大片看看| 亚洲真实伦在线观看| 亚洲av成人精品一二三区| 国产欧美亚洲国产| 日本午夜av视频| 青春草亚洲视频在线观看| 乱系列少妇在线播放| 老女人水多毛片| 日本-黄色视频高清免费观看| 女人精品久久久久毛片| 久久韩国三级中文字幕| 成年av动漫网址| .国产精品久久| 91久久精品国产一区二区三区| 99热这里只有精品一区| 久久99精品国语久久久| 亚洲欧美日韩另类电影网站| 三级国产精品欧美在线观看| 一区二区三区四区激情视频| 亚洲丝袜综合中文字幕| 亚洲一级一片aⅴ在线观看| 午夜福利影视在线免费观看| 国产亚洲一区二区精品| 亚洲婷婷狠狠爱综合网| 亚洲av综合色区一区| 久久久久久久久久久免费av| 美女国产视频在线观看| 老司机影院成人| 中文资源天堂在线| 亚洲欧美精品专区久久| 一区二区三区免费毛片| 欧美一级a爱片免费观看看| freevideosex欧美| 亚洲高清免费不卡视频| 亚洲av综合色区一区| 青春草视频在线免费观看| 久久久亚洲精品成人影院| 偷拍熟女少妇极品色| 国产精品一区www在线观看| 国产精品久久久久久精品古装| 99国产精品免费福利视频| 免费观看a级毛片全部| 嫩草影院入口| 久久狼人影院| 国产亚洲5aaaaa淫片| 久久影院123| 久久99热6这里只有精品| 久久这里有精品视频免费| 我的老师免费观看完整版| 2018国产大陆天天弄谢| 少妇熟女欧美另类| 岛国毛片在线播放| 日韩精品有码人妻一区| 亚洲av免费高清在线观看| 尾随美女入室| 亚洲国产精品专区欧美| 国产淫片久久久久久久久| 色5月婷婷丁香| 看非洲黑人一级黄片| 在现免费观看毛片| 国产成人午夜福利电影在线观看| 视频区图区小说| 久久97久久精品| 涩涩av久久男人的天堂| 久久热精品热| 日本av手机在线免费观看| 久久久久久人妻| 人体艺术视频欧美日本| 大码成人一级视频| 亚洲国产精品成人久久小说| 精品久久久噜噜| 亚洲精品aⅴ在线观看| 国产亚洲午夜精品一区二区久久| 国产亚洲精品久久久com| 一级爰片在线观看| 欧美日韩视频精品一区| av国产精品久久久久影院| 国产日韩欧美视频二区| 国产成人午夜福利电影在线观看| 久久热精品热| 26uuu在线亚洲综合色| 男人和女人高潮做爰伦理| 97超碰精品成人国产| 久久久久精品久久久久真实原创| 高清欧美精品videossex| 日日摸夜夜添夜夜添av毛片| 亚洲伊人久久精品综合| 曰老女人黄片| 亚洲国产日韩一区二区| 少妇的逼好多水| 婷婷色av中文字幕| 精品国产露脸久久av麻豆| 九色成人免费人妻av| 免费看av在线观看网站| 成人午夜精彩视频在线观看| 欧美精品人与动牲交sv欧美| 有码 亚洲区| 香蕉精品网在线| 欧美日韩一区二区视频在线观看视频在线| 亚洲自偷自拍三级| 欧美精品亚洲一区二区| 日韩成人伦理影院| 成年美女黄网站色视频大全免费 | av专区在线播放| 国产永久视频网站| 国产欧美另类精品又又久久亚洲欧美| 日本欧美国产在线视频| 亚洲第一av免费看| 国产永久视频网站| 久久人人爽人人爽人人片va| 午夜福利网站1000一区二区三区| 这个男人来自地球电影免费观看 | 久久久久久久久久久丰满| 亚洲无线观看免费| 日本与韩国留学比较| 人人妻人人添人人爽欧美一区卜| 国产一级毛片在线| 精品少妇黑人巨大在线播放| 亚洲欧洲国产日韩| 一区二区三区乱码不卡18| 成人美女网站在线观看视频| 男女免费视频国产| 菩萨蛮人人尽说江南好唐韦庄| 色视频在线一区二区三区| 九九爱精品视频在线观看| 岛国毛片在线播放| 我的女老师完整版在线观看| 日本欧美国产在线视频| 午夜免费观看性视频| 最新的欧美精品一区二区| 久久99热这里只频精品6学生| 国产精品一二三区在线看| 久久av网站| 免费播放大片免费观看视频在线观看| 能在线免费看毛片的网站| 观看免费一级毛片| 亚洲av国产av综合av卡| 美女福利国产在线| 乱码一卡2卡4卡精品| 纵有疾风起免费观看全集完整版| 免费av中文字幕在线| 亚洲国产色片| 97在线人人人人妻| 国产一区二区三区综合在线观看 | 国产片特级美女逼逼视频| 国产午夜精品一二区理论片| 免费人妻精品一区二区三区视频| 久久ye,这里只有精品| 久久久久久久大尺度免费视频| 日日啪夜夜撸| 国产精品熟女久久久久浪| 日韩,欧美,国产一区二区三区| 多毛熟女@视频| 久久久久久久久久久丰满| 欧美精品一区二区大全| 国产亚洲午夜精品一区二区久久| 亚洲人与动物交配视频| 亚洲欧美日韩东京热| 99热国产这里只有精品6| 最近中文字幕高清免费大全6| 建设人人有责人人尽责人人享有的| 欧美激情极品国产一区二区三区 | 99久国产av精品国产电影| 久热这里只有精品99| 汤姆久久久久久久影院中文字幕| 91在线精品国自产拍蜜月| 久久久久久久国产电影| 亚洲精品一二三| 国产极品天堂在线| 亚洲情色 制服丝袜| 精品亚洲乱码少妇综合久久| 天堂中文最新版在线下载| 自拍偷自拍亚洲精品老妇| 最近手机中文字幕大全| 国产深夜福利视频在线观看| 22中文网久久字幕| 我要看日韩黄色一级片| 亚洲国产成人一精品久久久| 性高湖久久久久久久久免费观看| 精品国产乱码久久久久久小说| 六月丁香七月| 青春草视频在线免费观看| 汤姆久久久久久久影院中文字幕| 91在线精品国自产拍蜜月| 丝袜喷水一区| 国产精品无大码| 嫩草影院入口| 成人综合一区亚洲| 深夜a级毛片| 六月丁香七月| 人妻 亚洲 视频| 亚洲美女黄色视频免费看| 日本av手机在线免费观看| 亚洲成人一二三区av| 99re6热这里在线精品视频| 视频区图区小说| 久久久久国产精品人妻一区二区| 人体艺术视频欧美日本| 亚洲国产精品一区三区| 日日啪夜夜爽| 免费高清在线观看视频在线观看| 久久精品久久精品一区二区三区| 日韩人妻高清精品专区| 亚洲欧洲精品一区二区精品久久久 | 国产亚洲欧美精品永久| 亚洲国产色片| 精品少妇黑人巨大在线播放| 桃花免费在线播放| 狂野欧美激情性bbbbbb| 亚洲性久久影院| 亚洲精品国产色婷婷电影| 蜜桃久久精品国产亚洲av| 大陆偷拍与自拍| 成人影院久久| 老熟女久久久| 少妇人妻一区二区三区视频| 亚洲国产色片| 日日摸夜夜添夜夜添av毛片| 91aial.com中文字幕在线观看| 亚洲欧美中文字幕日韩二区| 大香蕉久久网| 色吧在线观看| 久久久久久久久久人人人人人人| a级片在线免费高清观看视频| 国产一区亚洲一区在线观看| 国产乱来视频区| 熟女电影av网| 精品少妇久久久久久888优播| 国产淫语在线视频| 3wmmmm亚洲av在线观看| 曰老女人黄片| 免费大片黄手机在线观看| 啦啦啦在线观看免费高清www| 免费看日本二区| 午夜福利网站1000一区二区三区| 国产免费视频播放在线视频| 国产免费又黄又爽又色| 婷婷色综合大香蕉| 在线 av 中文字幕| 赤兔流量卡办理| 视频区图区小说| 精品一区二区三区视频在线| 日韩熟女老妇一区二区性免费视频| 国产亚洲午夜精品一区二区久久| 新久久久久国产一级毛片| 日本与韩国留学比较| 99热6这里只有精品| 在线看a的网站| 少妇裸体淫交视频免费看高清| 中文字幕人妻熟人妻熟丝袜美| 免费人妻精品一区二区三区视频| 久久青草综合色| kizo精华| 精品午夜福利在线看| 18禁动态无遮挡网站| 国产女主播在线喷水免费视频网站| 久久 成人 亚洲| 80岁老熟妇乱子伦牲交| 久久久久久伊人网av| 晚上一个人看的免费电影| 亚洲av国产av综合av卡| 少妇猛男粗大的猛烈进出视频| 97超碰精品成人国产| 爱豆传媒免费全集在线观看| 久久久久久久国产电影| 亚洲精品日韩在线中文字幕| 亚洲综合精品二区| 欧美激情国产日韩精品一区| 国产欧美另类精品又又久久亚洲欧美| 热re99久久精品国产66热6| 交换朋友夫妻互换小说| 亚洲欧美日韩卡通动漫| 另类亚洲欧美激情| 日韩大片免费观看网站| 晚上一个人看的免费电影| 国内少妇人妻偷人精品xxx网站| 日日爽夜夜爽网站| 久热这里只有精品99| 亚洲综合色惰| 成人二区视频| 在线看a的网站| 国产成人a∨麻豆精品| 欧美老熟妇乱子伦牲交| 亚洲一级一片aⅴ在线观看| 亚洲精品中文字幕在线视频 | 欧美精品人与动牲交sv欧美| 久久久久久久精品精品| 亚洲人成网站在线观看播放| 国产精品人妻久久久影院| 三上悠亚av全集在线观看 | 王馨瑶露胸无遮挡在线观看| 女性生殖器流出的白浆| 亚洲精品国产成人久久av| 乱人伦中国视频| 国产熟女午夜一区二区三区 | 熟女人妻精品中文字幕| 中文字幕人妻丝袜制服| 国产成人免费无遮挡视频| 亚洲精品一二三| 多毛熟女@视频| 成年女人在线观看亚洲视频| 菩萨蛮人人尽说江南好唐韦庄| 纯流量卡能插随身wifi吗| 久久99蜜桃精品久久| 乱系列少妇在线播放| 亚洲经典国产精华液单| 久久99热这里只频精品6学生| 国产精品不卡视频一区二区| 一级毛片电影观看| 久久久a久久爽久久v久久| 日韩不卡一区二区三区视频在线| 国产亚洲一区二区精品| 国产永久视频网站| 国产精品熟女久久久久浪| 国产av码专区亚洲av| 精品久久久精品久久久| 免费av不卡在线播放| 国产成人免费无遮挡视频| 精品酒店卫生间| 久久久久久久精品精品| 好男人视频免费观看在线| 日韩熟女老妇一区二区性免费视频| 日本爱情动作片www.在线观看| 午夜免费男女啪啪视频观看| 日日撸夜夜添| av线在线观看网站| 少妇人妻久久综合中文| 中文字幕免费在线视频6| 少妇人妻一区二区三区视频| 免费在线观看成人毛片| freevideosex欧美| 国产精品久久久久久精品电影小说| 人体艺术视频欧美日本| 国模一区二区三区四区视频| 免费看光身美女| 久久影院123| 18+在线观看网站| 欧美另类一区| 丁香六月天网| 最近最新中文字幕免费大全7| 高清在线视频一区二区三区| 91aial.com中文字幕在线观看| 最近的中文字幕免费完整| 精品熟女少妇av免费看| av国产久精品久网站免费入址| 三级经典国产精品| 国产av国产精品国产| 秋霞伦理黄片| 又粗又硬又长又爽又黄的视频| 老女人水多毛片| 美女国产视频在线观看| 亚洲精品一区蜜桃| 免费人妻精品一区二区三区视频| 一级黄片播放器| 一区二区三区精品91| 国产无遮挡羞羞视频在线观看| 最新的欧美精品一区二区| 99热这里只有是精品在线观看| 女人久久www免费人成看片| 毛片一级片免费看久久久久| 国产男女超爽视频在线观看| 精品久久久精品久久久|