The study of land use change is essential due to its significant effects on the environment and human life. Land use change modelers have mostly focused on the binary methods (e.g. urban and non-urban) rather than multiple land use changes methods. Also, most of the models used for modeling of land use changes are global parametric models (e.g. artificial neural network) and local non-parametric models (e.g. Multivariate Adaptive Regression Spline (MARS)) is rarely used to simulate multiple LUCs.  Local models split the data into subsets and fit distinct models on each of the subsets. Non-parametric models do not have a fixed model structure or model structure is unknown before the modeling. On the other hand, global models perform modeling using all the available data. In addition, parametric models have a fixed structure before the modeling. In this paper, we applied one of the well-known data mining tools, called multivariate adaptive regression spline, as one of the local non-parametric models with geospatial information system and satellite images to simulate urban and agriculture land use changes for northern part of Iran including cities of Sari and Ghaem Shahr over a period of 22 years during 1992 and 2014. Landsat images are the core source for information extraction and modeling of land use change in this research. Landsat images of 1992 (TM) and 2014 (ETM+) were used for modeling the urban and agricultural land uses changes. The spatial predictors considered for urban and agriculture modeling in this area were distance to urban areas, distance to agriculture areas, distance to roads, distance to water, aspect, and slope in 1992. After the modeling, a sensitivity analysis was performed on the effective parameters of the land use changes. The results of the sensitivity analysis verified that the most important factors were distances from agricultural and urban areas as well as elevation. To assess the model performance, the receiver operating characteristics (ROC) and total operating characteristics (TOC) were used. Considering multiple thresholds, ROC reveals how strong each threshold of the generated index is in diagnosing either presence or absence of a characteristic which results in a two by two contingency table without informing the size of each entries. While preserving the important information revealing by ROC, TOC gives size information of each entry. The area under the receiver operating characteristics curve for urban and agricultural land uses were 65% and 61.01%, respectively. Also, we have labeled thresholds for 0.67 and 0.40 in total operating characteristic curves for agriculture and urban gain to show four entries in the two-by-two contingency tables, respectively. These thresholds represent the probability of land use change for pixels in the suitability maps. According to the results, the percent of observations that are reference change and have been diagnosed as change by the model are equals to 36.8% and 67.06% for these thresholds, respectively.