Predicting the Potential Distribution of Hylomecon japonica in China under Current and Future Climate Change Based on Maxent Model

: Hylomecon japonica is considered a natural medicinal plant with anti ‐ inflammatory, anti ‐ cancer and antibacterial activity. The assessment of climate change impact on its habitat suitability is important for the wild cultivation and standardized planting of H. japonica . In this study, the maximum entropy model (Maxent) and geographic information system (ArcGIS) were applied to predict the current and future distribution of H. japonica species, and the contributions of variables were evaluated by using the jackknife test. The area under the receiver operating characteristic curve (AUC) value confirmed the accuracy of the model prediction based on 102 occurrence records. The predicted potential distributions of H. japonica were mainly concentrated in Jilin, Liaoning, Shaanxi, Chongqing, Henan, Heilongjiang and other provinces (adaptability index > 0.6). The jackknife ex ‐ periment showed that the precipitation of driest month (40.5%), mean annual temperature (12.4%), the precipitation of wettest quarter (11.6%) and the subclass of soil (9.7%) were the most important factors affecting the potential distribution of H. japonica . In the future, only under the shared socio ‐ economic Pathway 245 (SSP 245) scenario model in 2061–2080, the suitable habitat area for H. japon ‐ ica is expected to show a significant upward trend. The area under other scenarios may not increase or decrease significantly.


Introduction
Hylomecon japonica (Thunb.) Prantl et Kundig, a medicinal plant of the Papaveraceae family, is widely distributed in Shaanxi, Liaoning, Chongqing, Shanxi, Sichuan, Heilongjiang, and Henan in China, along with Japan and Korea [1]. It is a kind of perennial herb with irregular serrated leaves, and characterized by bright yellow flowers [2]. Meanwhile, its roots have medicinal potential, so it is commonly used in folk to treat diseases such as arthritis, neuralgia and eczema [3]. Various alkaloids, phenols, flavonoids, saponins and other active compounds isolated from H. japonica have been reported. A new species of endophytic bacteria belonging to sphingomonas was isolated from the rhizome of H. japonica [4][5][6]. Although many pharmacological studies, including anti-inflammatory, anticancer and antibacterial activities, have been reported [7][8][9], the investigation of the medicinal value and mechanism of H. japonica is far from enough. Due to the important medicinal value of the species, the wild resources of H. japonica have been greatly destroyed, and in combination with climate change, the natural habitat of these plants may gradually diminish. Most of the previous studies focused on the chemical composition of H. japonica, while other studies focused on the mechanism of its pharmacological action. However, up to now, little information is available on the guidance to protect, develop and utilize the H. japonica resources.
The studies for regionalization of traditional Chinese medicine (TCM) have been performed to analyze the relationship between the distribution of TCM resources and their surrounding environmental variables. Relevant factors such as ecological environment, geographical distribution, regional characteristics and quality of the TCM were classified. The research on regionalization of TCM began in the 1990s. With the rapid development of science and technology, the research methods of TCM regionalization were gradually improved [10]. The Species Distribution Model (SDM) is the main model used to predict suitable habitats for species in the regionalization of TCM. At present, there are a large number of SDMs and related software that can draw suitable habitat distribution maps of species based on the observed distribution data. The potential distribution of Nelumbo nucifera Gaertn. was assessed in China by the genetic algorithm for rule set production (GARP) and Maxent [11]; based on CLIMEX ecological software, the data were compiled and processed to predict the adaptive area of Veratrum taliense Loes. F [12]. Additionally, the BIOCLIM niche model has been used to predict the potential adaptability of Yulania Liliflora (Desr.) D. L, Fu [13]. Among these models, the Maxent model has shown many advantages over other models when applied to "existence-only" species occurrence data, and is widely used for species distribution prediction and habitat suitability assessment [14][15][16][17].
Combined with some other software such as ArcGIS, Maxent can extract factor variables according to the geographic information of sample collection points to analyze the habitat suitability of species. Furthermore, potentially suitable habitats can be divided into a specified number of levels as needed. Many reports have showed that the Maxent model has been used to analyze the habitat suitability of medicinal plants and ecologically important species in the study of regionalization of TCM. Based on the Maxent model and GIS technology, the production regionalization of Angelica in China was studied [18]. By constructing a Maxent model, Wang Dan predicted that the high adaptability areas of Bupleurum marginatum were mainly distributed in the border of provinces in the southern region [19]. The spatial analysis function of ArcGIS and Maxent model were used to predict the ecological and quality suitable areas of Gentianae Macrophyllae Radix in China [20]. The Maxent and ArcGIS were also used to predict the production regionalization and its suitability level of Paris polyphylla Smith var. chinensis (Franch.) Hara. in China [21]. The Maxent model only needs a set of known data (such as longitude and latitude information of sampling points) and variable factors, such as topography, precipitation, soil, temperature, vegetation type. These continuous or classified data are used to predict the habitat suitability distribution by combining the interactions between variables. The model has the characteristics of good performance with incomplete datasets, high efficiency, easy operation, high accuracy, small sample size requirements and high guiding significance [22].
Seeking suitable habitats and evaluating the main ecological factors that affect the distribution of H. japonica has a great practical significance for the rational development, utilization and protection of the resources. In the present study, the occurrence records of H. japonica were collected for habitat suitability assessment. We used Maxent to simulate the suitability of H. japonica habitats on the occurrence records, based on bioclimatic variable factors, and to predict the future distribution of H. japonica according to SSP 245 and SSP 585 climate change scenarios by means of ArcGIS. The Maxent model evaluated the accuracy of model calculation results according to receiver-operating characteristic (ROC) curve (AUC), and the key environmental variables related to the geographical distribution of H. japonica were mapped by jackknife test.

Species Data
In this study, a total of 102 occurrences of H. japonica were collected from the following sources: (1) the Global Biodiversity Information Facility Data Portal [23]; (2) the Specimen Resources Sharing Platform for Education [24]; (3) the Chinese Virtual Herbarium [25]; (4) Chinese plant species information System [26]; (5) relevant literature [2,6,27,28] and field investigation. In addition, only one was retained when the obtained geographic information of the sample repeated. When the obtained records lacked the required geographic information, Google Earth 7.0 was used to obtain approximate latitude and longitude information based on the described geographic location. The geographic information of each sampling point was converted into geographic coordinates (World Geodetic System 1984 data) through ArcGIS (version 10.2) [29]. Finally, there were a total of 102 valid points, and the sample points were mainly from Liaoning, Shaanxi, Chongqing, Shanxi, Sichuan and Heilongjiang in China. According to the requirements of Maxent software, the distribution records of the H. japonica species were converted into ʺ.csvʺ format files. The vector map of China's administrative divisions with province boundaries provided by the National Basic Geographic Information System [30] was used as the base map for analysis, and the species' suitability distribution map was generated by ArcMap 10.2.

Variables
Through relevant literature [31][32][33], 54 ecological factor data, including meteorological factors, soil type and topographic data, which affect the distribution of Chinese medicinal materials were selected (Table 1). Among them, meteorological factors played a key role in influencing the distribution of target plants, including 19 bioclimatic variables, precipitation variables and temperature variables, with a spatial resolution of 30 s. All of them, were downloaded from the WorldClim-Global Climate Database (WorldClim version 2.1) [34]. The BCC-CSM2-MR (the second-generation National (Beijing) Climate Center moderate resolution climate system model) of the sixth phase of the Coupled Model Intercomparison Projects (CMIP6) was selected as the future climate model of China. In the climate prediction experiment, the shared socioeconomic paths (SSPs) and representative concentration paths (RCP) were considered, and two scenarios SSP 2-RCP 4.5 (SSP 245) and SSP 5-RCP 8.5 (SSP 585) were selected for simulation. SSP 245 is the updated scenario of RCP 4.5, which uses a moderate development path, represents an intermediate level of greenhouse gas emissions and a nominal 4.5 Wm −2 radioactive forcing level by 2100 [35]; while, SSP 585 is the updated scenario of RCP 8.5, which uses the development path dominated by fossil fuels, represents a high level of greenhouse gas emissions and a nominal 8.5 Wm −2 radioactive forcing level by 2100 [36]. Additionally, soil factors have also been shown to be important factors affecting species distribution [37]. Thus, the soil variables which were obtained from Harmonized World Soil Database [38] were used to establish the model. Finally, the Digital Terrain Model (DTM) with 1 km resolution, downloaded from the data center of resources and environment science of Chinese Academy of Sciences [39], was used to generate the slope, aspect and elevation data layers (Table 1). In the projection coordinate system as "Asia_Lambert_Conformal_Conic", these variables were extracted into China and converted to ASCII files according to the unified processing scope and pixel size.

Classification Factors Name Unit
Meteorological factors bio1 Mean annual temperature °C bio2 Mean of monthly (max temperature-min temperature) °C bio3 Isothermality (bio2/bio7x100) 1 bio4 Standard deviation (SD) of temperature seasonality 1 bio5 Max temperature of warmest month °C bio6 Min temperature of coldest month °C bio7 Temperature annual range (bio5-bio6) °C bio8 Mean temperature of wettest quarter °C bio9 Mean temperature of driest quarter °C bio10 Mean temperature of warmest quarter °C bio11 Mean temperature of coldest quarter °C bio12 Annual precipitation mm bio13 Precipitation of wettest month mm bio14 Precipitation of driest month mm bio15 Coefficient of variation of precipitation seasonality 1 bio16 Precipitation of wettest quarter mm bio17 Precipitation of driest quarter mm bio18 Precipitation of warmest quarter mm bio19 Precipitation of coldest quarter mm Precipitation factors pre1-pre12 Monthly precipitation from January to December mm The temperature factors tmean1-tmean12 Average monthly temperature from January to December °C

Maximum Entropy (Maxent) Model
Since it was written in JAVA language by Phillips et al. [40] based on the principle of maximum entropy in 2006, Maxent software has been widely used to model species distribution [19]. The software implements the maximum entropy algorithm to generate the probability distribution of pixels in the modeled area grid. The core idea of maximum entropy theory is that the state of things with maximum entropy is the most likely to reflect or be close to the real state of things [41]. Maxent is a multivariable method, which selects the distribution with the maximum entropy from the qualified distribution as the optimal distribution, constrained by incomplete information about the distribution [42]. That is, Maxent needs a complete record set of the spatial distribution of species to construct the model, and then finds the constraint conditions (environmental variables) that limit the distribution of species, and then establishes the relationship between them to obtain the best estimate of species distribution and calculate the relative contribution rates of different environmental variables [43].
The method of the Maxent model to deal with the problem is to transform the problem being explored into a probability model, so that the randomness of the problem is expressed as a probability distribution, and the solution of the problem is to find the optimal probability distribution [41]. The calculation result of entropy value increases with the input of environmental factors related to each distribution data and the increase of iteration times. Finally, the state with maximum entropy is obtained, that is, the state closest to the reality of things [44]. Mathematically, given a random variable ε, it has n different potential results X1, X2... Xn, the occurrence probabilities are p1, p2... pn, respectively, the entropy of ε can be expressed by the following formula [45]: 1 Therefore, the distribution model of the species can be generated by inputting a set of environmental variables and a set of geospatial occurrence locations, and then Maxent can generate the probability map of species existence in three output formats: primitive, cumulative and logical [46]. We generally use a logical output format that obtains the estimated probability of existence for a given location, and the suitability value ranges from 0 (unsuitability) to 1 (optimum suitability) [46]. To further distinguish the suitability, thresholds can be used to reclassify the probability map, converting continuously predicted raster pixel values to binary values (0/1) and eventually displaying them as percentage data.
The advantage of the Maxent model is that it only needs the geographical spatial distribution data of the current distribution of the target species when simulating and predicting the distribution of the species. Moreover, it still has a good prediction effect when the correlation between species distribution points and environmental factors is not clear and the species distribution records are few [47]. Therefore, it is very suitable for species distribution modeling and habitat suitability prediction.

Modeling Procedure
Since there are many significant correlations between the extracted environment variables, the multicollinearity of variables leading to model overfitting should be avoided. Therefore, in establishing models, the contribution rate of each variable was evaluated through the knife-cut method of the Maxent model and the variables with small contributions were eliminated. Then, the correlation analysis of the remaining data ( Table 2) was performed using SPSS software (version 18.0, IBM Corp., Armonk, USA). When the correlation coefficient of two variables was less than 0.80, all the relevant variables were retained, and when the correlation coefficient of two variables was greater than 0.80, the variables with more important ecological significance were retained. For instance, the variables bio1 and bio9 were correlated (r = 0.931), considering the relative contribution result, bio9 was dropped. Finally, a total of 12 variables were retained for further analyses (Table 3).

Accuracy Test
The area under receiver-operating characteristic (ROC) curve (AUC) obtained by the accuracy test of the ROC curve analysis method was used in the Maxent model performance evaluation. The AUC values were between 0 and 1 and divided into six parts. When the AUC value was lower than 0.5, the model executed was worse than contingency. When the AUC value ranges from 0.5 to 0.6, the model performance was poor, 0.6-0.7 was fair, 0.7-0.8 was good, 0.8-0.9 was very good and 0.9-1 was excellent [50]. Figure  1 shows that the AUC value predicted by Maxent based on the variables of potential H. japonica distribution was 0.937, which indicates that Maxent model has an excellent prediction effect on the potential distribution region of H. japonica in China.

Contribution of Environmental Variables
The relative contributions of variables to the Maxent model under current circumstances obtained by Maxent iterative calculation and the results of jackknife analyses obtained the main variable factors influencing the H. japonica distribution and their contribution rates (Figure 2). The obtained variables affecting the distribution of H. japonica were as follows: precipitation of driest month (bio14), mean annual temperature (bio1), precipitation of wettest quarter (bio16), subclass of soil (sym90), texture class name and code (usda) and mean temperature of wettest quarter (bio8). The total contribution of these variables was 86.8%, indicating that the H. japonica distribution was strongly influenced by these variables. The total contribution rate of the bioclimatic variables was 72.8%. The total contribution rate of the soil variables was 16.8%. The contribution rate of the three terrain variables was 10.4%. This indicates that the influence of bioclimatic variables on H. japonica habitat distribution is much greater than that of topography and soil. The relationship between the precipitation in the driest month ( Figure 3A) and the probability of the existence of H. japonica indicated that when the precipitation in the driest month was about 5 mm, the probability of the existence of H. japonica could reach 80%. When precipitation in the driest months was between 4 and 35 mm, the probability of H. japonica was higher than 50%. The mean annual temperature was also an important bioclimatic variable which affects the distribution of H. japonica ( Figure 3B). When the mean annual temperature ranged from −7.5 to 7.5 °C, the presence probability of H. japonica occurrence was greater than 50%, when the mean annual temperature was in the range of −2-4 °C, the presence probability of H. japonica could reach 80% or more. Then, precipitation of wettest quarter also exerted an enormous function on the distribution of H. japonica ( Figure 3C). When the precipitation of wettest quarter ranged between 375 and 640 mm, the presence probability of H. japonica occurrence was greater than 50%, when the precipitation of wettest quarter was 400-490 mm, the presence probability of H. japonica could reach 80% or more. Finally, mean temperature of wettest quarter was another bioclimatic variable that affected the distribution of H. japonica ( Figure 3F). When the mean temperature of wettest quarter ranged from 19 to 27, more than 50% probability of H. japonica could be achieved. Among all of the soil variables, subclass of soil played a vital role in the distribution of H. japonica. In the Acrisols, Glossic Chernozems, Luvic Calcisols, Petric Calcisols, Gleyic Greyzems and Haplic Gypsisols ( Figure 3D) soils the presence probability of H. japonica were more than 50%. Topsoil USDA texture classification also played an important role in influencing the distribution of H. japonica. In the silt loam, loam and loamy sand ( Figure  3E) soils the presence probability of H. japonica was more than 50%. Figure 4 shows the habitat suitability distributions of H. japonica in China, according to the combination technology of Maxent and ArcGIS software. The habitat suitability results were expressed as probability with a range of 0-1 [51]. Using the reclassification tool of ArcMap 10.2, the probability results obtained were divided into four levels, of which 0-0.2 was considered unsuitable, 0.2-0.4 was considered generally suitable, 0.4-0.6 was considered moderately suitable and 0.6-1 was considered as highly suitable [52]. The analysis results showed that the regions highly suitable for H. japonica were mainly distributed in the southeast of Heilongjiang, Jilin, Liaoning and Shandong Provinces, as well as the southeast of Shaanxi, Shanxi, Gansu. Additionally, the junction of Shaanxi with Henan, Sichuan, Hubei and Chongqing, with relatively small distributions in Guizhou, Hunan, Jiangxi and Fujian, and the area would be predicted to be 7.03×10 5 km 2 , accounting for 35.22% of the total suitable area. The moderately suitable areas for H. japonica were found to be mainly distributed around the highly suitable areas of H. japonica, and the area would be predicted to be 4.92×10 5 km 2 , accounting for 24.64% of the total suitable area. In addition, the generally suitable areas for H. japonica were found to be distributed in Guizhou, Sichuan, and Heilongjiang Province, with sporadic distributions in Anhui, Jiangsu, Zhejiang, Henan, and a small part of Xinjiang and Xizang (Figure 4), and the total habitat area would be predicted to be 19.96×10 5 km 2 .

Future Changes in Suitable Habitat Areas of the H. japonica Species
The predicted future suitable climatic distributions of H. japonica under the SSP 245 and SSP 585 climate change scenarios for the periods of 2041-2060 and 2061-2080 are shown in Figure 5. Our results showed that during 2041-2060, for SSP 245, the highly suitable habitat area of H. japonica would be predicted to be 7.12×10 5 km 2 , the moderately suitable habitat area would be predicted to be 4.77×10 5 km 2 and the total habitat area would be predicted to be 19.88×10 5 km 2 . Under SSP 585, the highly suitable habitat area for H. japonica would be predicted to be 7.02×10 5 km 2 , the moderately suitable habitat area would be predicted to be 5.00×10 5 km 2 and the total would be predicted to be 19.81×10 5 km 2 . The main distribution regions and area are not much different from the current potential distribution ( Figure 5). During 2061-2080, under SSP 245, the highly suitable habitat area for H. japonica would be predicted to be 8.29×10 5 km 2 , mainly located in the southeast of Heilongjiang, Jilin, Liaoning and Shandong provinces, and the south of Shaanxi, Shanxi, as well as the east of Gansu, Yunnan and most parts of Chongqing; the moderately suitable habitat area would be predicted to be 5.21×10 5 km 2 and the total habitat area would be predicted to be 21.20×10 5 km 2 . Under SSP 585, the highly suitable habitat area for H. japonica would be predicted to be 7.13×10 5 km 2 , the moderately suitable habitat area would be predicted to be 4.93×10 5 km 2 and the total habitat area would be predicted to be 20.00×10 5 km 2 . In the future, only under the SSP 245 scenario model in 2061-2080, the suitable habitat area for H. japonica is expected to show a significant upward trend. The area under other scenarios will not increase or decrease significantly, that is, there is basically no change compared with the current one ( Figure 5).

Discussion
Traditional Chinese medicines, as traditional medicinal plants in China, play a pivotal role in the countryʹs medical and health services. The scientific research of Chinese medicinal materials is different from the research of general chemical and biological medicine. Its survival rate and quality are the result of multiple variables, which are related to many ecological factors (such as temperature, rainfall, soil, topography). Therefore, blind cultivation and planting when the suitability of the production area is not clear may cause problems such as a decrease in the survival rate of the medicinal materials and a decrease in the main medicinal components. In recent years, the study for regionalization of TCM has become more and more important in the study of the distribution and cultivation of Chinese medicine. The Maxent model combined with ArcGIS has developed into an important tool for evaluating the influence of various variables on species distribution in recent years. They are based on species existence records and various variable factors that may affect the species distribution and analyze their correlations to predict the current potential distribution of species. The regionalization studies of various medicinal plants such as Alsophila spinulosain [47], Camptotheca acuminata [14], Ammopiptanthus [53] and Rhodiola kirilowii [54] have been carried out successively, and have shown high accuracy. H. japonica is a valuable medicinal plant with a wide range of pharmacological effects such as antibacterial, anti-inflammatory, anticancer and other pharmacological functions, and its material basis is mainly a variety of alkaloids. This is the first analysis and simulation of the distribution of H. japonica based on the Maxent model, which can provide a basis for the protection and cultivation of H. japonica.
In present study, 12 variables combined with 102 existence records were systematically used to establish the Maxent models of H. japonica species. The results of the prediction model showed that the habitat suitability of H. japonica agreed strongly with the data set of H. japonica occurrence [55]. The jackknife test of variables revealed that precipitation of driest month, mean annual temperature, precipitation of wettest quarter, subclass of soil, Topsoil USDA texture classification and mean temperature of wettest quarter had high contributions to the distribution of H. japonica. The precipitation of driest month showed the highest contribution (40.5%), followed by the mean annual temperature with the contribution of 12.4%. The precipitation of wettest quarter, which contributed 11.6%, was another important variable to explain the distribution of H. japonica, and subclass of soil, which contributed 9.7%, was another important variable. However, under the SSP 245 and SSP 585 climate change scenarios, total annual precipitation is expected to increase slightly at a rate of 6.4 and 8.0 mm per decade, according to studies conducted by Li et al. [35]. Therefore, according to our experimental results, part of the reason why the suitable distribution area of H. japonica does not change significantly in the future under the two scenarios can be explained. Finally, from the results of the study, the characteristics of the ecological environment factors in the best area of H. japonica distribution are the precipitation of driest month (5 mm), the mean annual temperature (−2-4 °C), the precipitation of wettest quarter (400-490 mm) and the subclass of soil (Glossic Chernozems, Gleyic Greyzems, Haplic Gypsisols). These results are consistent with the habitat characteristics of the existing distribution of H. japonica [56]. The results indicated that H. japonica species could be introduced into many undiscovered potential regions, such as Shandong and Guizhou.
From the predictions, it was found that the highly suitable habitat for H. japonica increased during 2061-2080, under the SSP 245 climate scenario, from 7.65% to 9.03% in China. Under each climatic scenario, the highly suitable area of H. japonica will continue to be mainly concentrated in the current potential distribution area, and the overall pattern was consistent with the past distribution. For example, under the most changing climate scenario (2061-2080, SSP 245), the total suitable area of H. japonica has only increased by 6.30% compared to the current one, so it is predicted that the overall pattern will only change slightly. In addition, in the SSP 245 (2061-2080) climate scenario, the suitable distribution area of H. japonica extended from Qinling Mountains and Changbai Mountains to surrounding areas. This is consistent with global observations that climate warming promotes the migration of plants to high altitude and high latitude regions [57]. Researchers believe that global warming may have different effects on the distribution of different species by causing them to expand, shift and even contract [58,59]. However, this effect is uncertain. First, human activities and climate change may also promote the adaptation of plant species to new climatic conditions. Moreover, with the different and continuous changes in environmental protection policies in various countries, the trend, extent and speed of global warming are also uncertain. Therefore, analyzing the possible changes in the distribution of species allows us to develop strategies earlier to prevent or reverse the negative trends that may occur in the future.
This study can provide a scientific basis for the cultivation and planting of H. japonica in China. However, before large-scale planting, field trials must be carried out because the yield and quality of medicinal plants can be affected by many other factors, such as existing species in the introduced area, natural disasters and road traffic. According to our results, this method could be used to predict the potential distribution of other medicinal plants and provide a valuable tool for species conservation and distribution research. However, there are still some aspects to be improved in this study. First of all, the potential distribution area predicted by the model may overestimate the actual niche of the species, because the species may be unable to disperse to the possible areas due to human interference, topographical obstacles and species competition. Secondly, in addition to the selected variables in this study, other variables may also affect the distribution range of H. japonica. Additionally, in the future distribution prediction, this study includes soil and topographic data, although the accuracy of the results is increased on the premise that the future topography and soil environment will not change. Therefore, limitations are also increased. Finally, the combination of multi-effect components with quality evaluation will make the study of H. japonica regionalization more comprehensive and reasonable.

Conclusions
In this study, a large number of occurrence data sets for H. japonica species were obtained from China, and a potential distribution model of H. japonica species based on multiple variables was successfully established. The current and future habitat suitability of H. japonica under SSP 245 and SSP 585 scenarios in 2041-2060 and 2061-2080 was successfully simulated, and the main variables affecting its distribution were determined. The results indicated that precipitation of driest month, mean annual temperature, precipitation of wettest quarter, subclass of soil, Topsoil USDA texture classification and mean temperature of wettest quarter were important factors which could control the distribution of H. japonica species. This provides a reference for the protection, introduction and cultivation of H. japonica in ecologically suitable areas. However, in addition to the variables that were included, other factors (including biological factors, vegetation types and human activities) may also affect the survival rate and the main medicinal ingredient of medicinal plants, and may be key factors in the formation of metabolites of medicinal plants. Further studies involving more factors in the model will enhance the accuracy of the prediction of H. japonica distribution in China.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.