The bond interface between fiber-reinforced polymer (FRP) and masonry units is one of the weakest links in tensile strengthened masonry structures. Several empirical formulae have been proposed for estimation of bond strength between FRP reinforcements and masonry units. However, the accuracy of existing formulae for predicting bond strength seems to be significantly limited. In the present study, the M5´ and Multivariate adaptive regression splines approaches are employed to predict the bond strength between FRP reinforcement and masonry units. To develop new models, a comprehensive database including 575 test series (230 distinctive specimens) is collected from different sources in the literature. The newly proposed formulations consider several preeminent parameters involved in the debonding process including the reinforcement width, the ratio between widths of FRP reinforcement and masonry unit, the tensile strength of substrate, the axial strength of reinforcement, and bond length. A comparative study is conducted to evaluate the performances of the developed models against the well-known equations. Results indicated that the proposed models remarkably outperform the existing equations in terms of accuracy. Furthermore, a sensitivity analysis is done to determine the most important parameters in predicting the bond strength. Finally, the safety of different methods is evaluated based on the demerit point classification scale.