It is shown that the proton formation probabilities, extracted from experimental decay half-lives, can be well reproduced by a simple multilinear formula with only three parameters. The parameters obtained by considering the standard root mean square deviation and the mini-max criteria are very similar to each other. In addition, we applied Bayesian analysis to study the uncertainties of the parameters and the model predictions. In this way we explain the systematics of proton decay half-lives. The multilinearity of the model also provides a way to classify the relative hindrance of different proton decays. All the recent experimental data agree very well with the model prediction. Our Bayesian analysis suggests that those new data do help constrain the uncertainty of the model parameters.