The momentum space subtraction (MOM) scheme is one of the most frequently used renormalization schemes in perturbative QCD (pQCD) theory. In the paper, we make a detailed discussion on the gauge dependence of the pQCD prediction under the MOM scheme. Conventionally, there is renormalization scale ambiguity for the fixed-order pQCD predictions, which assigns an arbitrary range and an arbitrary error for the fixed-order pQCD prediction. The principle of maximum conformality (PMC) adopts the renormalization group equation to determine the magnitude of the coupling constant and hence determines an effective momentum flow of the process, which is independent to the choice of renormalization scale. There is thus no renormalization scale ambiguity in PMC predictions. To concentrate our attention on the MOM gauge dependence, we first apply the PMC to deal with the pQCD series. We adopt the Higgs boson decay width, $Gamma(Hto gg)$, up to five-loop QCD contributions as an example to show how the gauge dependence behaves before and after applying the PMC. It is found that the Higgs decay width $Gamma (Hto gg)$ depends very weakly on the choices of the MOM schemes, being consistent with the renormalization group invariance. It is found that the gauge dependence of $Gamma(Hto gg)$ under the $rm{MOMgg}$ scheme is less than $pm1%$, which is the smallest gauge dependence among all the mentioned MOM schemes.