In this study using Frequency Ratio (FR), Statistical Index (SI), Weights Of Evidence(WOF), Logistic Regression (LR), Random Forest (RF) models the probability of gully formation was calculated in Aytamar watershed and susceptibility maps was prepared. First the thematic maps of 13 gully conditioning factors including lithological formations, distance to faults, faults density, altitude, slope-length, slope angle, slope aspect, plan curvature, profile curvature, distance to roads, land use, distance to rivers, stream power index and topographic wetness index was prepared. Then landslide inventory map was combined with each gully conditioning factor and all models weights and parameters were calculated. Area under curve for test data was calculated as 0.74, 0.78, 0.75, 0.86 and 0.96 for Frequency Ratio (FR), Statistical Index (SI), Weights Of Evidence(WOF), Logistic Regression (LR), Random Forest (RF) models, respectively. Random forest, Frequency Ratio and Logistic Regression have the least the area of high susceptibility zone, respectively. With respect three validation criteria multivariate methods including Random Forest and Logistic Regression had the best performance among all models.