Iran is the second largest country in the west of Asia, with a population of approx. 80 million people. This county comprises 31 provinces and has a total area of 1.648 million square kilometers, of which 8.7% is cultivated land. Average annual precipitation is 228 mm and approx. 90% of the country is arid or semi-arid. In 2014, main crops included wheat, barley, rice and maize, which were grown on 53%, 16%, 5% and 2% of the country’s cultivated areas, respectively. For empirical study in this work, Iran’s 31 provinces were considered. In each province, cropping pattern of four dominant crops (wheat, barley, rice, and maize) were simulated and analysed for the next 25 years. Most of the farms in Iran are monoculture and in-rotation cultivated areas are negligible, then in this study we assumed that all agricultural areas will be cultivated monoculture.

A large amount of data, including national and provincial data, was used in this research ( Table 1). Three types of data have been used: (1) some data were input directly into the SD sub-model for demand estimation; (2) some data were used by the crop allocation sub-model to simulate the distribution of estimated demand within the provinces; and (3) some data were used only for calibration of the model.

The integrated model consists of two sub-models: SD sub-model and crop allocation sub-model. SD estimates the demanded area for agricultural land at national level by considering the interaction between agronomics driving forces and without spatial resolution. Crop allocation sub-model simulates distribution of demand within geographical regions of the country by considering spatial characteristics of each region. Fig.1 shows the overall structure of the model. The exogenous variables of the model include population growth, economic growth, climate change, and ecological features of each region. The output of the model is presented in two steps. First at country scale, the model estimates the demand, the production, the trade and the area of farmland of each crop. Second at regional scale, the model determines the share of each crop from the total cultivated area of each region.

Sterman (2000) developed five steps to create SD models: 1) Problem articulation; 2) Dynamic hypothesis; 3) Formulation and deployment in software; 4) Testing; 5) Policy formulation and evaluation.

Recently, Mesgari et al . (2017) have presented a SD model to predict the total agricultural land use demand at the national scale. In this study, a modified version of that model has been used to estimate the demand, production, trade, and cultivated area of each crop at national scale. Fig. 2 describes the causal diagram of this SD model. The major assumptions in the development of this sub-model include: (1) demand for agricultural products will increase as population growth and income growth; (2) increasing demand will stimulate the supply side through the price mechanism; (3) part of the demand will be satisfied through the domestic production and the rest is met by trade mechanisms; (4) to achieve equilibrium, policy makers may apply optimization policies to interact between demand and supply such as: import tariffs, export incentives, and subsidies.

The stock-flow diagram (SFD) is the core of SD model, and is the process of quantization and materialization of the causal loop diagram ( Li et al. , 2012 ). Fig. 3 shows the SFD of the model which has been depicted in the Vensim software. Three main feedback loops could be distinguished in this diagram:

(1) Demand loop: In this loop, increasing " GDP " and increasing " population " lead to higher " demand " for crops. GDP and population are exogenous variables and are the main drivers of demand . Consequently, increasing demand cause an increase in the " price " of crops. Then the increased price will reduce demand and balance it in a negative loop.

(2) Supply loop: In this loop, rising price stimulate and increase " production " in a positive loop due to rising " value added per worker " and enhancing " cultivated land area " consequently. Increasing production will reduce price and smooth it.

(3) Trade loop: In this loop, a drop in price (due to decreased demand or increased production ) will reduce " ratio of in/out price " and then, import will be economically justified. Increasing import will reduce price and smooth it. On the contrary, this mechanism would be applied to export . In fact, trade in agricultural products creates equilibrium between prices within the country and abroad.

The main output of SD sub-model is the demanded area for each crop at country level. This output is the input of crop allocation sub-model. Complete details about the variables and formulas are available in Table S1 [suppl].

By inspiration from Dyna-CLUE ( Verburg & Overmars, 2009 ), an algorithm has been suggested to project the share of each region from national demanded cultivated area for each crop (referred in this paper as "crop allocation"). Dyna-CLUE allocates a land use to each cell based on competition between different land uses while crop allocation procedure determines portion of each region area which will be dedicate to a specific crop cultivation based on competition between different regions. In this algorithm, as depicted in Fig. 4, first we must calculate the total available land for each crop in each region based on total agriculturally usable area, spatial policies, allowable conversions and other specific restrictions. Then we must calculate the relative capability or advantage of each region for the cultivation of each crop. Finally, algorithm distributes the demanded area that is derived from SD sub-model among regions taking into account the total available land and the relative capabilities of each region.

More specific, the crop allocation sub-model like any other allocation mechanism needs to consider the following three parts: (1) cultivated land requirements (demand); (2) available land area (supply) by considering the spatial policies and restrictions; (3) allocation mechanism based on region capabilities.

(1) Demand: the area demanded to cultivate different crops will be determined by the SD sub-model described in the previous section.

(2) Supply: the total available land in each region should be calculated based on spatial policies, restrictions, conversion rules and area of protected lands. For this, we must reduce the total area of restricted lands, protected lands, and non-convertible lands from the total area of agriculturally usable lands.

(3) Allocation mechanism based on capabilities: it is assumed that crop allocation depends on the suitability of each region and current cropping pattern. In fact, crops tend to be cultivated in more suitable areas in long term due to rationality of farmers and ecological enforcements. But these changes could not be happen suddenly because of farmers’ adhesion to current situation and their resistance to change and the need for time to develop infrastructure. Then two factors influence on the future cropping pattern: (1) suitability and (2) change elasticity. In this study, two criteria have been introduced to represent these two factors (see Eq. 1 and Eq. 2) inspired by Balassa index ( Balassa, 1965 ).

The data of 1996 to 2006 were used to calibrate the SD sub-model in terms of parameters specifying and variable settings. Econometric techniques were used to obtain the response coefficients of SD model ( e.g. price response of demand, productivity response to production). In cases where estimation techniques could not be applied ( e.g. due to short time series), response coefficients have been embedded from the literature. Initial value of stock variables in SD model will be equal to its value in the first year of simulation.

The main parameter that needs to be calibrated in the crop allocation sub-model is suitability coefficient of each crop. By taking 1998 as the base year, value of parameter C was determined so that the most correlation occurs between actual cropping pattern and simulated cropping pattern in 2006. Results show that the values of C in Iran and for an 8 years’ time step is 0.3, 0.4, 0.5 and 0.3 for wheat, barley, rice and maize, respectively.

There are two ways to validate a SD model: structural validation and behavioral validation ( Barlas, 1996 ). Structural validation is the process of determining whether the theories and assumptions underlying the conceptual model are correct, reasonable, and sufficient. Behavioral validation is the numerical process of determining whether the model behaviors accurately represent the conceptual framework of the model ( Qudrat-Ullah, 2012 ; Mesgari et al. , 2017 ). Since we developed the SD model based on the expert knowledge and the available literature, structural validation of the model has been approved in several focus groups in academic centers related to agricultural and economic studies. Also the model’s behavioral validation is demonstrated by comparing its results to historical trends using data of 2007 to 2014. According to Barlas (1996) , a model will be valid if the standard error rate is smaller than 5%. In this research, examined variables include production and cultivated area . Table 2 reports standard errors for the period 2007-2014. Also Fig. 5 shows the dynamic trend of historical values of cultivated area and production variables, compared with their simulated values for wheat and barley between 2007-2014. All the error rates were smaller than 5%, which confirms the validity of the proposed SD model.

In addition, the data from2006 to 2014 was used to verify the whole integrated model and to evaluate its ability for projection by comparing its results to historic trends. R-squared values of the regression model between the historical data and simulated data were 0.80, 0.76, 0.61, and 0.76 ( p <0.01) for wheat, barley, rice and maize, respectively. These results indicate that the model is reliable and it could be used to project the crop pattern in the future based.

In this study, one probable scenario is defined for the future 25 years in accordance with the combinations of exogenous variables including population growth (medium), economic growth (4.9 %) and precipitation amount (212 mm). This scenario is based on the world population prospects of the United Nations and it almost keeps the present pace with demographic and economic development. Other scenarios could be designed and compared which we leave them to future studies.

Empirical results could be presented in two categories. First, the total cultivated area of each crop that is the output of SD sub model. Second, the spatial distribution of cultivated areas in the geospatial regions. Fig. 6 presents the changes in the area of national cultivated land for wheat, barley, rice, and maize during 2015–2040. It is obvious that although a constant growth in total area of cultivated land is predicted, there is a significant difference in the intensity of variation for different crops. There are almost 25 million hectares of agriculturally usable land in Iran, of which only 16 million hectares are currently under cultivation. Then the predicted increase for cultivated land of different crops will be mainly from two sources: agricultural development in agricultural lands that are not currently under cultivation, and displacement of crop’s farms.

The total area of cultivated land represents a considerable and stable growth for wheat and barley, where the area is predicted to be about 7.7 and 2.2 million hectares in 2040 for wheat and barley, respectively. Most of this increase occurs during the first period of simulation of the model and then, it will continue very slowly. Wheat shows the highest growth rates among the four crops, while the total area of rice and maize will be approximately constant.

Fig. 7 illustrates the national geospatial distribution of cultivated areas for four crops in 2015 and 2040 and in 31 provinces. These maps represent the share of each crop from the total agricultural land area of the province, i.e. , the range of 25-50% for wheat means that in the province, 25-50% of the total agricultural land is allocated to wheat. The characteristics of the predicted changes in geospatial distribution of major crops could be interpreted from Fig. 6 and Fig. 7.

Wheat is the most important crop in Iran, traditionally cultivated throughout the country. The growth rate of wheat cultivated lands will be higher in the western mountainous regions and less in the central arid regions. Just like wheat, barley is also a geographically ubiquitous crop; its cultivation area will be the second largest in Iran. Rice is mainly grown in northern Iran, where high temperature and high precipitation are favor to rice cultivation. Of course, some rice can be distinguished in other regions due to historic and the desire of farmers to maintain the previous state. Maize farms have the lowest share of cereal cultivated areas, and the area of its lands will decrease in all provinces. This may be due to the water shortage in Iran in the future.