In this paper, an implementation of Backtracking search optimization (BSA), a non-gradient iterative evolutionary algorithm, is presented for model parameter estimation applications using gravity anomalies. This is the first study to use this promising bio-inspired and population-based global optimization algorithm to invert geophysical data sets. We concentrated on demonstrating the efficiency and robustness of the proposed algorithm using some gravity anomalies originated from anticlinal and synclinal structures. Anomaly equations of these geological structures are based on planar approximation. We restricted our study to relatively simple profile data sets. The nature of the inverse problem was assessed by producing cost function topography image maps for the model parameters. The resolvability characteristics and unanticipated possible dependencies between the model parameter pairs were determined using these low misfit region maps. Parameter tuning studies enabled us to determine optimum values for the algorithm-based control parameters of the BSA, and to get the best efficiency from the algorithm. Some synthetic anomalies and two field gravity data sets observed over the Pays De Bray anticline system (France) were successfully inverted through BSA algorithm. In addition, a Markov Chain Monte Carlo method performing Metropolis-Hastings sampler revealed that all solutions obtained are within the confidence intervals without uncertainties. This study shows that the BSA algorithm can be an alternative inversion tool to the most widely applied global optimizers for geophysical anomalies, such as other evolutionary-based algorithms and particle swarm optimization. Additionally, BSA can be used for more complex multi-dimensional geophysical problems having formidable anomaly equations.